maginit.c

Go to the documentation of this file.

/*
    This file is part of magpar.

    Copyright (C) 2002-2009 Werner Scholz

    www:   http://www.magpar.net/
    email: magpar(at)magpar.net

    magpar is free software; you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation; either version 2 of the License, or
    (at your option) any later version.

    magpar is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with magpar; if not, write to the Free Software
    Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
*/

static char const Id[] = "$Id: maginit.c 2681 2009-07-31 04:30:53Z scholz $\n\n";
static char const Td[] = "$Today: " __FILE__ "  " __DATE__ "  "  __TIME__  " $\n\n";

#include "init.h"

int MagInit(GridData *gdata)
{
  MagparFunctionLogBegin;

  int rank,size;
  ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
  ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);

  /* TODO: rely on readkrn to broadcast n_prop -
     we need it as along as M is a shared vector during serial init.
  */

  PetscPrintf(PETSC_COMM_WORLD,"Creating Vec M for %i nodes\n",gdata->n_vert);
  PetscSynchronizedPrintf(PETSC_COMM_WORLD,
    "<%i>local nodes: %i  n_prop: %i\n",
    rank,gdata->ln_vert,gdata->n_prop
  );
  PetscSynchronizedFlush(PETSC_COMM_WORLD);

  ierr = VecCreate(PETSC_COMM_WORLD,&gdata->M);CHKERRQ(ierr);
  ierr = VecSetSizes(gdata->M,ND*gdata->ln_vert,PETSC_DECIDE);CHKERRQ(ierr);
  ierr = VecSetFromOptions(gdata->M);CHKERRQ(ierr);
  ierr = VecSetBlockSize(gdata->M,ND);CHKERRQ(ierr);

  int       *init_mag;
  PetscReal *init_magparm;
  ierr = PetscMalloc(gdata->n_prop*sizeof(int),&init_mag);CHKERRQ(ierr);
  ierr = PetscMalloc(gdata->n_prop*sizeof(PetscReal),&init_magparm);CHKERRQ(ierr);

  for (int i=0;i<gdata->n_prop;i++){
    init_mag[i]=3;
    init_magparm[i]=0.0;
  }

  int nmax = gdata->n_prop;

  PetscTruth flg;
  ierr = PetscOptionsGetIntArray(PETSC_NULL,"-init_mag",(PetscInt*)init_mag,(PetscInt*)&nmax,&flg);CHKERRQ(ierr);
  if (flg!=PETSC_TRUE) {
    init_mag[0]=3;
    PetscPrintf(PETSC_COMM_WORLD,
      "Option -init_mag not found, using default value: %i\n",
      init_mag[0]
    );
  }
  if (nmax<=1) {
    ierr = PetscOptionsGetReal(PETSC_NULL,"-init_magparm",init_magparm,&flg);CHKERRQ(ierr);
    if (flg!=PETSC_TRUE) {
      init_magparm[0]=0;
      PetscPrintf(PETSC_COMM_WORLD,
        "Option -init_magparm not found, using default value: %g\n",
        init_magparm[0]
      );
    }

    ierr = MagSet(-1,init_mag[0],init_magparm[0],gdata);CHKERRQ(ierr);
  }
  else if (nmax==gdata->n_prop) {
    nmax = gdata->n_prop;
    ierr = PetscOptionsGetRealArray(PETSC_NULL,"-init_magparm",init_magparm,(PetscInt*)&nmax,&flg);CHKERRQ(ierr);
    if (flg!=PETSC_TRUE) {
      PetscPrintf(PETSC_COMM_WORLD,
        "Option -init_magparm not found, using default value %g for all volumes\n",
        init_magparm[0]
      );
    }
    else if (nmax!=gdata->n_prop) {
      SETERRQ2(PETSC_ERR_ARG_INCOMP,
        "Wrong number of values for option -init_magparm! %i != %i\n",
        nmax,gdata->n_prop
      );
    }
    for (int i=0;i<gdata->n_prop;i++){
      ierr = MagSet(i,init_mag[i],init_magparm[i],gdata);CHKERRQ(ierr);
    }
  }
  else {
    SETERRQ2(PETSC_ERR_ARG_INCOMP,
      "Wrong number of values for option -init_mag! %i != %i\n",
      nmax,gdata->n_prop
    );
  }

  ierr = PetscFree(init_mag);CHKERRQ(ierr);
  ierr = PetscFree(init_magparm);CHKERRQ(ierr);

  MagparFunctionLogReturn(0);
}