00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024 static char const Id[] = "$Id: main.c 2962 2010-02-04 19:50:44Z scholz $\n\n";
00025 static char const Td[] = "$Today: " __FILE__ " " __DATE__ " " __TIME__ " $\n\n";
00026
00027 #ifdef PYTHON
00028 #include "Python.h"
00029 #endif
00030
00031 #include "init/init.h"
00032 #include "util/util.h"
00033 #include "io/magpario.h"
00034 #include "field/field.h"
00035
00036 #ifdef SUNDIALS_VERSION
00037 #include "llg/llg.h"
00038 #endif
00039
00040 #ifdef TAO
00041 #include "emini/emini.h"
00042 #endif
00043
00044 #ifdef ADDONS
00045 #include "addons/addons.h"
00046 #include "llgts/llgts.h"
00047 #include "llg2/llg2.h"
00048 #include "brown/brown.h"
00049 #include "util/util.h"
00050 #endif
00051
00052 #ifdef EBM
00053 #include "ebm/ebm.h"
00054 #endif
00055
00056
00057 int Solve(GridData *gdata);
00058
00059 int main(int argc,char **args)
00060 {
00061 int ierr;
00062 GridData gdata;
00063
00064 ierr = PetscInitialize(
00065 &argc,&args,
00066 "allopt.txt",
00067 "magpar - parallel micromagnetics package\n"
00068 );CHKERRQ(ierr);
00069 #ifdef TAO
00070 TaoInitialize(&argc,&args,(char *)0,"magpar - parallel micromagnetics package\n");
00071 #endif
00072
00073 #ifdef PYTHON
00074 Py_Initialize();
00075
00076
00077 #define NEWPYTHONIMPORT
00078 #ifdef NEWPYTHONIMPORT
00079 gdata.pModule = PyImport_ImportModule(MAGPARSTDPYTHONMOD);
00080 #else
00081 PyObject *pName;
00082 pName = PyString_FromString(MAGPARSTDPYTHONMOD);
00083 gdata.pModule = PyImport_Import(pName);
00084 #endif
00085 if(!gdata.pModule) {
00086 PetscPrintf(PETSC_COMM_WORLD,"Could not import Python module %s.py. The module file needs to be in PYTHONPATH!\n",MAGPARSTDPYTHONMOD);
00087 }
00088 else {
00089 PetscPrintf(PETSC_COMM_WORLD,"Imported Python module %s.py.\n",MAGPARSTDPYTHONMOD);
00090 }
00091 #endif
00092
00093 PetscLogDouble t_t1,t_t2;
00094 PetscPrintf(PETSC_COMM_WORLD,">> %s::%s\n",__FILE__,__FUNCT__);
00095 PetscGetTime(&t_t1);
00096
00097 PetscPrintf(PETSC_COMM_WORLD,
00098 "\nmagpar is distributed under the terms of the GNU General Public License\n"
00099 "Copyright (C) 2002-2010 Werner Scholz\n"
00100 );
00101
00102 #include "magpar_version.h"
00103 #include "magpar_revision.h"
00104 PetscPrintf(PETSC_COMM_WORLD,
00105 "--------------------------------------------------------------------------\n\n"
00106 "magpar version %s build %s\n",
00107 MAGPAR_VERSION,MAGPAR_REVISION
00108 );
00109
00110
00111
00112
00113
00114 ierr = PetscMallocDebug(PETSC_TRUE);CHKERRQ(ierr);
00115
00116 ierr = SerInit(&gdata);CHKERRQ(ierr);CHKMEMQ;
00117
00118 ierr = ParInit(&gdata);CHKERRQ(ierr);CHKMEMQ;
00119
00120
00121 ierr = PetscMallocDebug(PETSC_FALSE);CHKERRQ(ierr);
00122
00123 ierr = Solve(&gdata);CHKERRQ(ierr);
00124
00125 PetscGetTime(&t_t2);
00126 PetscPrintf(PETSC_COMM_WORLD,
00127 "<< %s::%s took %g s = %g h\n",
00128 __FILE__,__FUNCT__,t_t2-t_t1,(t_t2-t_t1)/3600
00129 );
00130
00131 #ifdef PYTHON
00132 #ifndef NEWPYTHONIMPORT
00133 Py_DECREF(pName);
00134 #endif
00135 if (gdata.pModule) {
00136 Py_DECREF(gdata.pModule);
00137 }
00138
00139 Py_Finalize();
00140 #endif
00141
00142 PetscFinalize();
00143 }
00144
00145
00146 int keepsolving(GridData *gdata)
00147 {
00148 MagparFunctionInfoBegin;
00149
00150
00151 static PetscReal jfinal=PETSC_MIN;
00152 static PetscReal mfinal[ND]={PETSC_MIN,PETSC_MIN,PETSC_MIN};
00153 static PetscReal tstop=PETSC_MIN;
00154
00155 #ifdef PETSC_USE_SINGLE
00156
00157
00158
00159 if (jfinal==(PetscReal)PETSC_MIN) {
00160 #else
00161 if (jfinal==PETSC_MIN) {
00162 #endif
00163 PetscTruth flg;
00164 ierr = PetscOptionsGetReal(PETSC_NULL,"-jfinal",&jfinal,&flg);CHKERRQ(ierr);
00165 if (flg!=PETSC_TRUE) {
00166 jfinal=-0.95;
00167 PetscPrintf(PETSC_COMM_WORLD,
00168 "Option -jfinal not found, using default value: %g\n",
00169 jfinal
00170 );
00171 }
00172 PetscPrintf(PETSC_COMM_WORLD,"jfinal: %g (1)\n",jfinal);
00173
00174 int nmax = ND;
00175 ierr = PetscOptionsGetRealArray(PETSC_NULL,"-mfinal",mfinal,(PetscInt*)&nmax,&flg);CHKERRQ(ierr);
00176
00177 if (flg!=PETSC_TRUE) {
00178 PetscPrintf(PETSC_COMM_WORLD,
00179 "Option -mfinal not found, using default value: %g,%g,%g\n",
00180 mfinal[0],mfinal[1],mfinal[2]
00181 );
00182 nmax=ND;
00183 }
00184 if (nmax!=ND){
00185 SETERRQ(PETSC_ERR_ARG_INCOMP,
00186 "Wrong number of values for option -mfinal!\n");
00187 }
00188 PetscPrintf(PETSC_COMM_WORLD,"mfinal: %g,%g,%g (1)\n",mfinal[0],mfinal[1],mfinal[2]);
00189
00190 ierr = PetscOptionsGetReal(PETSC_NULL,"-ts_max_time",&tstop,&flg);CHKERRQ(ierr);
00191 if (flg!=PETSC_TRUE) {
00192 if (gdata->mode==0 || gdata->mode==2 || gdata->mode==3) {
00193 tstop=5;
00194 }
00195 else {
00196 tstop=PETSC_MAX;
00197 }
00198 PetscPrintf(PETSC_COMM_WORLD,
00199 "Option -ts_max_time not found, using default value: %g\n",
00200 tstop
00201 );
00202 }
00203 tstop /= 1e9*gdata->tscale;
00204 PetscPrintf(PETSC_COMM_WORLD,"reduced end time: %g (1)\n",tstop);
00205 }
00206
00207
00208
00209
00210
00211
00212 int checkpoint=0;
00213
00214 static PetscLogDouble inp_last_CPUt=PETSC_MIN;
00215 ierr = PetscGetTime(&t_t2);CHKERRQ(ierr);
00216 if (t_t2-inp_last_CPUt > 2000.0) {
00217 checkpoint=1;
00218 }
00219
00220
00221
00222
00223 ierr = MPI_Bcast(&checkpoint,1,MPI_INT,0,PETSC_COMM_WORLD);CHKERRQ(ierr);
00224
00225
00226 if (checkpoint) {
00227 int oldinp;
00228
00229 oldinp=gdata->inp;
00230 gdata->inp=9999;
00231
00232
00233 PetscPrintf(PETSC_COMM_WORLD,
00234 "Writing checkpoint file (9999) at: %g ns\n",
00235 gdata->time*gdata->tscale*1e9
00236 );
00237 ierr = WriteSet(gdata);CHKERRQ(ierr);
00238
00239 gdata->inp=oldinp;
00240 inp_last_CPUt=t_t2;
00241
00242 checkpoint=0;
00243 }
00244
00245 int exitcode=1;
00246
00247 PetscReal valMpHext;
00248 ierr = MpHext(gdata->M,&valMpHext);CHKERRQ(ierr);
00249 if (valMpHext < jfinal) {
00250 PetscPrintf(PETSC_COMM_WORLD,
00251 "J//Hext %g <= jfinal %g\n",
00252 valMpHext,jfinal
00253 );
00254 exitcode=0;
00255 }
00256
00257 PetscReal Mtot[ND];
00258 ierr = Vec3VolAvg(gdata->M,Mtot);CHKERRQ(ierr);
00259
00260 for (int i=0;i<ND;i++) {
00261 if (Mtot[i] < mfinal[i]) {
00262 PetscPrintf(PETSC_COMM_WORLD,
00263 "M_%i %g <= mfinal\n",
00264 i,Mtot[i]
00265 );
00266 exitcode=0;
00267 }
00268 }
00269
00270 if (gdata->time > tstop) {
00271 PetscPrintf(PETSC_COMM_WORLD,
00272 "gdata->time %g ns >= tstop %g ns\n",
00273 gdata->time*gdata->tscale*1e9,tstop*gdata->tscale*1e9
00274 );
00275 exitcode=0;
00276 }
00277
00278
00279 ierr = MPI_Bcast(&gdata->mode,1,MPI_INT,0,PETSC_COMM_WORLD);CHKERRQ(ierr);
00280
00281 MagparFunctionInfoReturn(exitcode);
00282 }
00283
00284
00285 int Solve(GridData *gdata)
00286 {
00287 MagparFunctionLogBegin;
00288
00289 PetscPrintf(PETSC_COMM_WORLD,
00290 "--------------------------------------------------------------------------\n"
00291 "Entering solution loop\n\n"
00292 );
00293
00294 while (keepsolving(gdata) && gdata->mode!=99) {
00295 switch (gdata->mode) {
00296 case 0:
00297 #ifdef SUNDIALS_VERSION
00298 ierr = myTSStepPVode(gdata);CHKERRQ(ierr);
00299
00300
00301
00302
00303
00304
00305
00306
00307
00308
00309
00310 ierr = Htot_Energy(gdata);CHKERRQ(ierr);
00311 ierr = CheckIterationLLG(gdata);CHKERRQ(ierr);
00312 #endif
00313 break;
00314 case 1:
00315 #ifdef TAO
00316 ierr = EminiSolve(gdata);CHKERRQ(ierr);
00317 ierr = CheckIterationEmini(gdata);CHKERRQ(ierr);
00318 #if 0
00319
00320 gdata->mode=1;
00321 #endif
00322 #endif
00323 break;
00324 case 2:
00325 case 3:
00326 #ifdef SUNDIALS_VERSION
00327 ierr = myTSStepPVode(gdata);CHKERRQ(ierr);
00328 ierr = Htot_Energy(gdata);CHKERRQ(ierr);
00329 ierr = CheckIterationLLG(gdata);CHKERRQ(ierr);
00330 #endif
00331 break;
00332 #ifdef ADDONS
00333 case 4:
00334 ierr = myTSStep2(gdata);CHKERRQ(ierr);
00335 ierr = Htot_Energy(gdata);CHKERRQ(ierr);
00336 ierr = CheckIterationLLG2(gdata);CHKERRQ(ierr);
00337 break;
00338 #endif
00339 #ifdef EBM
00340 case 50:
00341 ierr = myTSStepPVodeEBM(gdata);CHKERRQ(ierr);
00342 ierr = CheckIterationEBM(gdata);CHKERRQ(ierr);
00343 gdata->equil=0;
00344 break;
00345 #endif
00346 #ifdef ADDONS
00347 case 51:
00348 ierr = myTSStep(gdata);CHKERRQ(ierr);
00349 ierr = Htot_Energy(gdata);CHKERRQ(ierr);
00350 ierr = CheckIterationLLGTS(gdata);CHKERRQ(ierr);
00351 break;
00352 case 52:
00353 ierr = BrownSolve(gdata);CHKERRQ(ierr);
00354 break;
00355 case 53:
00356 #ifdef TAO
00357 ierr = EminiSolve(gdata);CHKERRQ(ierr);
00358 ierr = CheckIterationEmini(gdata);CHKERRQ(ierr);
00359 #endif
00360 break;
00361 #endif
00362 default:
00363 SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Unknown mode selected! Exiting!\n");
00364 break;
00365 }
00366 }
00367
00368 PetscGetTime(&t_t2);
00369 PetscPrintf(PETSC_COMM_WORLD,
00370 "--------------------------------------------------------------------------\n"
00371 "Finished solution loop in %g s = %g h\n\n",
00372 t_t2-t_t1,(t_t2-t_t1)/3600
00373 );
00374
00375
00376
00377
00378
00379 gdata->inp=PetscAbsInt(gdata->inp);
00380 ierr = WriteLog(gdata);CHKERRQ(ierr);
00381 ierr = WriteSet(gdata);CHKERRQ(ierr);
00382
00383 MagparFunctionLogReturn(0);
00384 }
00385