Mathlink mit dem cernlib FORTRAN-Paket MINUIT verwenden

  • Ich konnte Mathlink erfolgreich zum Aufrufen einer einfachen Fortran-Subroutine verwenden. Als ich jedoch zu einem realistischen Problem kam, schlug ich fehl. Um genau zu sein, möchte ich einen Datensatz mit einer Funktion unter Verwendung des Cernlib-Pakets MINUIT (für Funktionsminimierung und Fehleranalyse) anpassen. Innerhalb von FORTRAN kann ich dann die besten Parametersätze und deren Fehler und viele mehr erhalten.

    Nun möchte ich einige Zahlen mit den resultierenden besten Parametern in Mathematica darstellen. Natürlich kann ich die Parameter mit den Augen ablesen ... Aber es wäre schön, dies mit Mathlink zu tun. Ich habe jedoch versagt. Ich habe eine ausführbare Datei erstellt, die von Mathlink jedoch nicht empfohlen werden kann. Installieren Sie einfach, bis ich es abgebrochen habe. Dies ist ein wenig merkwürdig, da der Fortran-Code in C erfolgreich aufgerufen wurde, als Mathlink nicht verwendet wurde. Hier ist der verkürzte FORTRAN-Code, der für ein MINUIT-Beispiel (noch etwas langwierig) myprog.f genannt werden soll = "http://www.slac.stanford.edu/comp/physics/cernlib/minuit/doc/minuit-at-SLAC.html" rel = "nofollow"> hier :

       ! n is not used, just to give the Mathlink function an argument
      subroutine myprog(n, fmin) 
      !program myprog  
      implicit none
      integer SYSIN/1/, SYSRD/5/, SYSWR/6/, SYSPU/7/
      integer n, ndat, npari, nparx, istat
      real(8) fmin, fedm, errdef
      external myfcn
    
      !     initialize unit numbers
      call mintio(SYSRD, SYSWR, SYSPU)
      !driven card running the fit
      open(unit=SYSRD, file='myprog.mincards', status='OLD') 
      !output, if commented out, the output will be written to terminal
      open(unit=SYSWR, file='myprog.out', status='UNKNOWN')      
    
      !     invoke minuit
      call minuit(myfcn, 0)      
      ! current status of minimization
      CALL mnstat(fmin, fedm, errdef, npari, nparx, istat)
      end
    
    
      subroutine myfcn(npar, g, f, x, iflag, futil)
      implicit none
      integer npar, iflag, i, ndat, n
      integer SYSIN/1/, SYSRD/5/, SYSWR/6/, SYSPU/7/
      parameter (ndat=5)
      real*8 f, g(*), x(*), xdat(ndat), ydat(ndat), edat(ndat), func,        
     &       answer/15.7/, rms, avg, diff, futil
      !     Save all data that must be preserved between calls!
      save xdat, ydat, edat         ! MUST
    
      if (iflag .eq. 1) then
      !     Initialization mode
         do i=1,ndat
            xdat(i) = i
            ydat(i) = answer + i - float(ndat + 1) / 2.
            edat(i) = 1.
         end do
      end if
      !     compute LSQ
      f = 0
      do i = 1, ndat
         f = f + (ydat(i) - func(xdat(i), x(1)))**2 / edat(i)
      end do
      end
    
      function func(x, par)
      implicit none
      real*8 x, par(*), func
      func = par(1)
      end
     

    Die getriebene Karte myprog.mincards , wo der Anfang steht Werte der Parameter werden gesetzt, enthält

     set title
    Minuit data cards for myprog example program
    parameters
    1 'average' 0.0 1.0
    
    set print 1
    mig
    mino
    return
    stop    
     

    Die Subroutine myprog(n,fmin) (fmin gibt den minimierten Wert von Chi-Quadrat an) kann in C mit folgendem Code erfolgreich aufgerufen werden:

     #include "stdio.h"
    #ifdef __cplusplus
    extern "C"{
    #endif
    
    double myprog(int n) 
    {
      double fmin;  
      myprog_(&n, &fmin);
      return fmin;
    }
    
    #ifdef __cplusplus
    }
    #endif
    
    int main(void)
    {
    printf("%f\n", myprog());
    return 0;
    }
     

    Das Ergebnis ist 10. Ich dachte, solange die FORTRAN-Subroutine in C aufgerufen werden kann, sollte man es in Mathematica mit Mathlink aufrufen können. Hier ist der C-Wrapper (ich füge den vollständigen Code ein, da er für andere hilfreich sein könnte):

     #include "mathlink.h"
    
    #ifdef __cplusplus
    extern "C"{
    #endif
    
    double myprog(int n) 
    {
      double fmin;
    
      myprog_(&n, &fmin);
      MLPutDouble(stdlink, fmin);
      return fmin;
    }
    
    #ifdef __cplusplus
    }
    #endif
    
    #if WINDOWS_MATHLINK
    
    #if __BORLANDC__
    #pragma argsused
    #endif
    
    int PASCAL WinMain( HINSTANCE hinstCurrent, HINSTANCE hinstPrevious, LPSTR lpszCmdLine, int nCmdShow)
    {
    char  buff[512];
    char FAR * buff_start = buff;
    char FAR * argv[32];
    char FAR * FAR * argv_end = argv + 32;
    
    hinstPrevious = hinstPrevious; /* suppress warning */
    
    if( !MLInitializeIcon( hinstCurrent, nCmdShow)) return 1;
    MLScanString( argv, &argv_end, &lpszCmdLine, &buff_start);
    return MLMain( (int)(argv_end - argv), argv);
    }
    
    #else
    
    int main(int argc, char* argv[])
    {
    return MLMain(argc, argv);
    }
    
    #endif
     

    Die Mathlink-Vorlagendatei myprog.tm lautet

     :Begin:
    :Function:       myprog
    :Pattern:        Minuit[n_Integer]
    :Arguments:      { n }
    :ArgumentTypes:  { Integer }
    :ReturnType:     Manual
    :End:
    
    :Evaluate: Minuit::usage = "Minuit[n] does fit using MINUIT. n can be an arbitary integer. This function returns the value of the least chi-square."
     

    Ich verwende Windows 7 mit MinGW. Ich habe math.bat geschrieben, um eine ausführbare Datei für Mathlink zu erstellen.

     gfortran -c myprog.f
    gcc -c callfmath.c -o myprogc.o
    mprep myprog.tm -o myprogtm.c
    gcc -c myprogtm.c -o myprogtm.o
    gcc *.o libminuit.a libML32i3.a -lm -lpthread -mwindows -lstdc++ -lgfortran -o myprog_math
    rm *.o
     
    19 May 2012
    unstable
1 answer
  • Hier beschreibe ich die Schritte, die ich unternommen habe, damit MacBook Pro unter Mac OS X 10.7.4 funktioniert. Unter Windows sollte es nicht sehr unterschiedlich sein, aber die Idee ist die gleiche. Ich werde versuchen, diesen Beitrag so einfach wie möglich zu gestalten und alle Details darzulegen, die ich kann.

    Fortran-Bibliotheken

    Jetzt beginnen indem Sie die Fortran-Bibliotheken kompilieren. Dazu müssen Sie die MINUIT -Routinen herunterladen. Ich werde an dem Verzeichnis /Users/jmlopez/Desktop/MathLink arbeiten. In diesem Verzeichnis habe ich ein Verzeichnis mit dem Namen fortran erstellt. Der MINUIT-Download wird in diesem Ordner fortran abgelegt.

    Wenn Sie den MINUIT-Ordner aus der ZIP-Datei extrahiert haben, werden viele Fortran-Dateien angezeigt. Eines wird insbesondere NEWd506cm.FOR genannt. Sie müssen den Namen dieser Datei in d506cm.F

    ändern. Navigieren Sie nun mit Ihrem Terminal (Mac und Linux) zu diesem Ordner und geben Sie den Befehl ein:

     gfortran -c *.F
     

    Damit werden nun alle Fortran-Dateien kompiliert. Der nächste Schritt besteht darin, die .o -Dateien in einer Bibliothek abzulegen. Dies können Sie durch Eingabe von

     ar rvs ../libminuit.a *.o
     

    tun. Dieser Befehl fügt alle .o -Dateien in der Datei libminuit zusammen. a, das sich in einem Verzeichnis über dem aktuellen Verzeichnis befindet. Jetzt benötigen wir Ihre Fortran-Datei, die diese Bibliothek verwendet. Ich werde es hier noch einmal posten, in der Hoffnung, dass, wenn jemand den Code kopiert und in eine Datei einfügt, es problemlos funktioniert, ich auch das Argument n loswird, da Sie es dort offensichtlich nicht brauchen. Hier ist die Datei myprog.f

           ! FILE: myprog.f
          subroutine myprog(fmin) 
          !program myprog  
          implicit none
          integer SYSIN/1/, SYSRD/5/, SYSWR/6/, SYSPU/7/
          integer ndat, npari, nparx, istat
          real(8) fmin, fedm, errdef
          external myfcn
    
          !     initialize unit numbers
          call mintio(SYSRD, SYSWR, SYSPU)
          !driven card running the fit
          open(unit=SYSRD, file='myprog.mincards', status='OLD')
          !output, if commented out, the output will be written to terminal
          open(unit=SYSWR, file='myprog.out', status='UNKNOWN')
    
          !     invoke minuit
          call minuit(myfcn, 0)
          ! current status of minimization
          CALL mnstat(fmin, fedm, errdef, npari, nparx, istat)
          end
    
          subroutine myfcn(npar, g, f, x, iflag, futil)
          implicit none
          integer npar, iflag, i, ndat, n
          integer SYSIN/1/, SYSRD/5/, SYSWR/6/, SYSPU/7/
          parameter (ndat=5)
          real*8 f, g(*), x(*), xdat(ndat), ydat(ndat), edat(ndat), func
          real*8 answer/15.7/, rms, avg, diff, futil
          !     Save all data that must be preserved between calls!
          save xdat, ydat, edat         ! MUST
    
          if (iflag .eq. 1) then
          !     Initialization mode
          do i=1,ndat
                xdat(i) = i
                ydat(i) = answer + i - float(ndat + 1) / 2.
                edat(i) = 1.
          end do
          end if
          !     compute LSQ
          f = 0
          do i = 1, ndat
                f = f + (ydat(i) - func(xdat(i), x(1)))**2 / edat(i)
          end do
          end
    
          function func(x, par)
          implicit none
          real*8 x, par(*), func
          func = par(1)
          end
     

    Jetzt können wir sie kompilieren Diese Datei, um das Objekt zu erhalten myprog.o

     gfortran -c myprog.f
     

    Wir sollten unsere Fortran-Bibliotheken verwenden. An dieser Stelle habe ich die Dateien myprog.o und libminuit.a im Verzeichnis /Users/jmlopez/Desktop/MathLink/fortran.

    C-Programm

    Um sicherzustellen, dass wir mit diesen Bibliotheken arbeiten können, erstellen Sie ein einfaches C pProgramm, das dies tut. Jetzt arbeiten wir im Verzeichnis /Users/jmlopez/Desktop/MathLink/cxx. In diesem Verzeichnis werde ich nun die Datei fortranCall.c erstellen:

     // FILE: fortranCall.c
    #include <stdio.h>
    
    double myprog_(double*);
    
    double myprog(void) {
      double fmin;  
      myprog_(&fmin);
      return fmin;
    }
    
    int main(void) {
        printf("%f\n", myprog());
        return 0;
    }
     
    < Beachten Sie, dass ich einen Prototyp für die Fortran-Funktion deklariert habe. Auf diese Weise wird der c-Compiler nicht beschweren, dass er nicht weiß, was myprog_ ist. Jetzt können wir diese Datei kompilieren, indem Sie Folgendes tun:

      gcc -lgfortran ../fortran/libminuit.a ../fortran/myprog.o fortranCall.c -o fortranCall
     

    Jetzt haben wir unsere Datei fortranCall. Wenn wir es ausführen, erhalten wir einen Fehler.

     Manuel-Lopezs-MacBook-Pro:cxx jmlopez$ ./fortranCall
    At line 13 of file myprog.f (unit = 5, file = 'stdin')
    Fortran runtime error: File 'myprog.mincards' does not exist
     

    I Beachten Sie hierbei, dass die MathLink-Verbindung mit Mathematica aus irgendeinem Grund abstürzt, weil Ihr Programm abstürzt. Als ich es zum ersten Mal tat, wurde die Verbindung hergestellt, aber sobald ich es lief, brach die Verbindung ab und ich hatte keine Ahnung, was los war, bis ich daran erinnerte, dass ich nicht die myprog.mincards -Datei hatte. Sobald wir diese haben Datei im selben Verzeichnis sollten wir folgendes erhalten:

     Manuel-Lopezs-MacBook-Pro:cxx jmlopez$ ./fortranCall
    10.000000
     

    Beachten Sie, dass das Programm die Datei generiert myprog.out. Am Ende der Datei steht:

     CALL TO USER FUNCTION WITH IFLAG = 3
    
    ..........MINUIT TERMINATED AND RETURNS TO USER PROGRAM. 
     

    Ich bin mir nicht sicher, ob dies gut oder schlecht ist Es bringt Probleme mit dem MathLink-Programm.

    Mathematica

    Ich werde jetzt in /Users/jmlopez/Desktop/MathLink/mma arbeiten. In diesem Ordner habe ich 5 Dateien: mlprog.m, MLFile.nb, mlprog.tm, myprog.mincards und mlprog.cpp.

    mlprog.cpp

     // FILE: mlprog.cpp
    #include "mathlink.h"
    
    extern "C" {
        void myprog_(double*);
    }
    
    void MyProg(void) {
        double fmin = 200.0; //myprog_(&fmin);
        MLPutReal64(stdlink, fmin);
    }
    
    int main(int argc, char* argv[]) {
        return MLMain(argc, argv);
    }
     

    Beachten Sie, wie ich den Teil kommentierte, in dem wir die Fortran-Routine aufrufen. Wir tun dies, um sicherzustellen, dass wir ein funktionierendes MathLink-Programm erhalten können. Diese Datei deklariert die Funktion MyProg. Diese Funktion akzeptiert keine Argumente und gibt nur einen Wert zurück: 200.

    mlprog.tm

     :Evaluate: BeginPackage["mlprog`"]
    
    :Evaluate: MyProg::usage = "MyProg[], returns a real number. "
    :Evaluate: MyProg::argerr = "This function takes no arguments.";
    :Evaluate: MyProg[args___] := (Message[MyProg::argerr]; $Failed)
    
    :Evaluate: EndPackage[]
    
    :Evaluate: Begin["mlprog`Private`"]
    
    :Begin:
    :Function: MyProg
    :Pattern: MyProg[]
    :Arguments: {}
    :ArgumentTypes: {Manual}
    :ReturnType: Manual
    :End:
    
    :Evaluate: End[]
     

    In dieser Vorlagendatei habe ich einige mma-Befehle deklariert, damit wir einige Eingabefehler beheben können. Da Sie keine Eingabe benötigen, schreiben Sie nichtsre.

    Schließlich benötigt die Datei mlprog.m nur eine Zeile:

     Install["mlprog"];
     

    Nun geht es weiter zum Terminal. Angenommen, wir befinden uns in dem Verzeichnis, in dem sich all diese Dateien befinden, mache ich nun Folgendes:

     export MMA=/Applications/Mathematica.app/SystemFiles/Links/MathLink/DeveloperKit/CompilerAdditions
     

    Ich mache das, da sich dieses Verzeichnis nicht in meinem Pfad befindet und ich derzeit zu faul bin, es dort aufzunehmen. Dieses Verzeichnis enthält die Binärdateien mprep und andere von mathematica in Mac OS X benötigte Dateien.

     ${MMA}/mprep mlprog.tm -o mlprogtm.cpp
    g++ -c -I${MMA} mlprogtm.cpp
    g++ -c -I${MMA} mlprog.cpp
    g++ -I${MMA} mlprogtm.o mlprog.o ../fortran/myprog.o ../fortran/libminuit.a -L${MMA} -lMLi3 -lstdc++ -lgfortran -framework Foundation -o mlprog
     

    Jetzt haben wir die ausführbare Datei mlprog. Hier ist ein Screenshot meines Mathematica-Notebooks, das diese ausführbare Datei verwendet:

    Mathlink mit dem cernlib FORTRAN-Paket MINUIT verwenden

    So weit so gut. Lassen Sie uns nun diesen Abschnitt in der Datei mlprog.cpp auskommentieren, um die fortran-Routine zu verwenden. Wir wiederholen alle obigen Schritte, um das Programm erneut zu kompilieren, und so läuft es für mich ab:

    Mathlink mit dem cernlib FORTRAN-Paket MINUIT verwenden

    Was ist hier gerade passiert? Die Funktion funktionierte beim ersten Mal, danach wurde die Verbindung abgebrochen. Wir haben gesehen, dass vor dem Einsatz der Fortran-Routine alles in Ordnung war. Und Sie haben erwähnt, dass Sie es mit einer anderen einfachen Fortran-Routine zum Laufen bringen konnten. Ich schätze, all das bedeutet, dass irgendwo in der Fortranroutine, hauptsächlich aus der Bibliothek, eine Anweisung zum Beenden des Programms vorhanden ist, die die Verbindung automatisch beendet. Dies sagt das Ende der erzeugten Datei aus:

      CALL TO USER FUNCTION WITH IFLAG = 3
    
     ..........MINUIT TERMINATED AND RETURNS TO USER PROGRAM.            
     ***************************************************************************
      MINUIT RELEASE 96.03  INITIALIZED.   DIMENSIONS 100/ 50  EPSMAC=  0.89E-15
     ***************************************************************************
                              MINUIT DATA BLOCK NO.   1
     ***************************************************************************
     stop                                              
     ******************************************************************************
     END OF DATA ON UNIT NO.  5
     TWO CONSECUTIVE EOFs ON PRIMARY INPUT FILE WILL TERMINATE EXECUTION.
    
     ..........MINUIT TERMINATED BY UNRECOVERABLE READ ERROR ON INPUT. 
     

    Es scheint, dass es beim ersten Aufruf gut funktioniert, aber beim zweiten Mal gibt es einen Fehler. Ich weiß nicht, wie ich das beheben kann, da ich die Fortran-Routinen nicht kenne, aber da haben Sie es. Ein Beispiel für die Verwendung von Fortran von Mathematica. Ich hoffe, dass Sie unter Windows keine großen Änderungen vornehmen werden.

    UPDATE

    I mal wieder im c-programm gesucht. Nachdem ich das Programm zweimal ausgeführt habe, bekomme ich dieselbe myprog.out-Datei. Nein

    18 May 2012
    Jon Galloway