parinit.c

Go to the documentation of this file.
00001 /*
00002     This file is part of magpar.
00003 
00004     Copyright (C) 2002-2010 Werner Scholz
00005 
00006     www:   http://www.magpar.net/
00007     email: magpar(at)magpar.net
00008 
00009     magpar is free software; you can redistribute it and/or modify
00010     it under the terms of the GNU General Public License as published by
00011     the Free Software Foundation; either version 2 of the License, or
00012     (at your option) any later version.
00013 
00014     magpar is distributed in the hope that it will be useful,
00015     but WITHOUT ANY WARRANTY; without even the implied warranty of
00016     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00017     GNU General Public License for more details.
00018 
00019     You should have received a copy of the GNU General Public License
00020     along with magpar; if not, write to the Free Software
00021     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00022 */
00023 
00024 static char const Id[] = "$Id: parinit.c 2962 2010-02-04 19:50:44Z scholz $\n\n";
00025 static char const Td[] = "$Today: " __FILE__ "  " __DATE__ "  "  __TIME__  " $\n\n";
00026 
00027 #include "init.h"
00028 #include "io/magpario.h"
00029 #include "field/field.h"
00030 
00031 #ifdef SUNDIALS_VERSION
00032 #include "llg/llg.h"
00033 #endif
00034 
00035 #ifdef TAO
00036 #include "emini/emini.h"
00037 #endif
00038 
00039 #ifdef PNG
00040 #include "png/writepng.h"
00041 #endif
00042 
00043 #ifdef ADDONS
00044 #include "addons/addons.h"
00045 #include "llgts/llgts.h"
00046 #include "llg2/llg2.h"
00047 #include "brown/brown.h"
00048 #endif
00049 
00050 #ifdef EBM
00051 #include "ebm/ebm.h"
00052 #endif
00053 
00054 int ParInit(GridData *gdata)
00055 {
00056   PetscTruth   flg;
00057 
00058   MagparFunctionLogBegin;
00059 
00060   PetscPrintf(MPI_COMM_WORLD,
00061     "--------------------------------------------------------------------------\n"
00062     "Starting parallel part\n\n"
00063   );
00064 
00065   /* initialize mode here: Hdemag needs it */
00066   ierr = PetscOptionsGetInt(PETSC_NULL,"-mode",(PetscInt*)&gdata->mode,&flg);CHKERRQ(ierr);
00067   if (flg!=PETSC_TRUE) {
00068     gdata->mode=99;
00069     PetscPrintf(PETSC_COMM_WORLD,
00070       "Option -mode not found, using default value: %i\n",
00071       gdata->mode
00072     );
00073   }
00074   PetscPrintf(PETSC_COMM_WORLD,"run mode: %i\n",gdata->mode);
00075 
00076   /* modify material properties */
00077   ierr = ModifyPropPar(gdata);CHKERRQ(ierr);
00078 
00079   ierr = EleVertVol(gdata);CHKERRQ(ierr);
00080 
00081   gdata->equil=0;
00082   /* Initialize with large value to avoid artificial equilibrium at start of simulation.
00083      Also done in first call to EquilCheck().
00084   */
00085   gdata->vequil=1.0;
00086 
00087   PetscReal  tstart;
00088   ierr = PetscOptionsGetReal(PETSC_NULL,"-ts_init_time",&tstart,&flg);CHKERRQ(ierr);
00089   if (flg!=PETSC_TRUE) {
00090     tstart=0;
00091     PetscPrintf(PETSC_COMM_WORLD,
00092       "Option -ts_init_time not found, using default value: %g\n",
00093       tstart
00094     );
00095   }
00096   PetscPrintf(PETSC_COMM_WORLD,"initial time: %g\n",tstart);
00097   tstart /= 1e9*gdata->tscale;
00098   PetscPrintf(PETSC_COMM_WORLD,"reduced initial time:      %g (1)\n",tstart);
00099   gdata->time=tstart;
00100 
00101   ierr = Htot(gdata);CHKERRQ(ierr);
00102   ierr = Htot_Energy(gdata);CHKERRQ(ierr);
00103 
00104   switch (gdata->mode) {
00105     case 0:
00106       /* LLG time integration */
00107 #ifdef SUNDIALS_VERSION
00108       ierr = myTSCreatePVode(gdata);CHKERRQ(ierr);
00109       ierr = WriteFEMAVS(gdata);CHKERRQ(ierr);
00110       ierr = CheckIterationLLG(gdata);CHKERRQ(ierr);
00111 #else
00112       SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"PVODE not compiled in! Exiting!\n");
00113 #endif
00114     break;
00115     case 1:
00116       /* energy minimization */
00117 #ifdef TAO
00118       ierr = myTSCreateEmini(gdata);CHKERRQ(ierr);
00119       ierr = WriteFEMAVS(gdata);CHKERRQ(ierr);
00120       ierr = CheckIterationEmini(gdata);CHKERRQ(ierr);
00121 #else
00122       SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"TAO not compiled in! Exiting!\n");
00123 #endif
00124     break;
00125     case 2:
00126       /* relax with high damping and no precession during t<0,
00127          then switch to normal damping
00128       */
00129 #ifdef SUNDIALS_VERSION
00130       ierr = myTSCreatePVode(gdata);CHKERRQ(ierr);
00131       ierr = WriteFEMAVS(gdata);CHKERRQ(ierr);
00132       ierr = CheckIterationLLG(gdata);CHKERRQ(ierr);
00133 #else
00134       SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"PVODE not compiled in! Exiting!\n");
00135 #endif
00136     break;
00137       /* relax with energy minimization, then switch to LLG time integration */
00138     case 3:
00139 #ifdef TAO
00140       ierr = myTSCreateEmini(gdata);CHKERRQ(ierr);
00141       ierr = WriteFEMAVS(gdata);CHKERRQ(ierr);
00142       ierr = CheckIterationEmini(gdata);CHKERRQ(ierr);
00143       PetscPrintf(PETSC_COMM_WORLD,">> EminiSolve: Relaxing magnetization using energy minimization\n");
00144       ierr = EminiSolve(gdata);CHKERRQ(ierr);
00145       gdata->inp=PetscAbsInt(gdata->inp);
00146 
00147 #ifdef SUNDIALS_VERSION
00148       ierr = myTSCreatePVode(gdata);CHKERRQ(ierr);
00149       ierr = CheckIterationLLG(gdata);CHKERRQ(ierr);
00150 #else
00151       SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"PVODE not compiled in! Exiting!\n");
00152 #endif
00153 
00154 #else
00155       SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"TAO not compiled in! Exiting!\n");
00156 #endif
00157     break;
00158 #ifdef ADDONS
00159     case 4:
00160       ierr = WriteFEMAVS(gdata);CHKERRQ(ierr);
00161       ierr = CheckIterationLLG2(gdata);CHKERRQ(ierr);
00162     break;
00163 #endif
00164     case 50:
00165 #ifdef EBM
00166       /* nudged elastic band method */
00167       ierr = myTSCreateEBM(gdata);CHKERRQ(ierr);
00168       ierr = WriteFEMAVS(gdata);CHKERRQ(ierr);
00169 #else
00170       SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"EBM not compiled in! Exiting!\n");
00171 #endif
00172     break;
00173 #ifdef ADDONS
00174     case 51:
00175       /* LLG time integration with PETSc internal ODE solvers */
00176       ierr = myTSCreate(gdata);CHKERRQ(ierr);
00177       ierr = WriteFEMAVS(gdata);CHKERRQ(ierr);
00178       ierr = CheckIterationLLGTS(gdata);CHKERRQ(ierr);
00179     break;
00180     case 52:
00181       ierr = BrownInit(gdata);CHKERRQ(ierr);
00182     break;
00183     case 53:
00184       /* energy minimization with non-linear magnet */
00185 #ifdef TAO
00186       /* shorten magnetization vectors to avoid starting in saturation */
00187       ierr = VecScale(gdata->M,0.1);CHKERRQ(ierr);
00188       /* need to recalculate field, energy after scaling M */
00189       ierr = Htot(gdata);CHKERRQ(ierr);
00190       ierr = Htot_Energy(gdata);CHKERRQ(ierr);
00191 
00192       ierr = myTSCreateEmini(gdata);CHKERRQ(ierr);
00193       ierr = WriteFEMAVS(gdata);CHKERRQ(ierr);
00194       ierr = CheckIterationEmini(gdata);CHKERRQ(ierr);
00195 #else
00196       SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"TAO not compiled in! Exiting!\n");
00197 #endif
00198     break;
00199 #endif
00200     case 99:
00201       /* just calculate fields, energies, etc. once and exit */
00202       ierr = WriteLog(gdata);CHKERRQ(ierr);
00203       ierr = WriteFEMAVS(gdata);CHKERRQ(ierr);
00204       ierr = WriteSet(gdata);CHKERRQ(ierr);
00205     break;
00206     default:
00207       SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"-mode out of valid range! Exiting!\n");
00208   }
00209 
00210   /* destroy various objects, which are not required any more */
00211   ierr = DataDestroyInit(gdata);CHKERRQ(ierr);
00212 
00213   MagparFunctionLogReturn(0);
00214 }

magpar - Parallel Finite Element Micromagnetics Package
Copyright (C) 2002-2010 Werner Scholz