Parallel RHS Vector

The RHS function still had to be implemented using PETSc code, and a number of computational tasks needed to be implemented, using the PETSc matrix-vector manipulation methods to calculate the various contributions to the total energy. The MicroMag code holds the magnetic field in two different data structures. The first is a two-dimensional array of size m(NMAX,3) where NMAX is the number of finite element vertices, and 3 represents the three physical dimensions in Cartesian space. The second data structure is a linearised one dimensional array, Y(3*NMAX)

The code switches between these two data structures as the Sundials solver expects a one dimensional array, but it is both more natural to code, to use the other data structure. In the main loop, the magnetic field is packed into the linearised array, and the solver called. The RHS function is called by the solver, and one of its arguments is Y . The RHS function then has to unpack this one-dimensional array back into the two dimensional array. The PETSc data structures called Vectors are not arrays in the FORTRAN sense, and can only be manipulated by calling the appropriate methods. However, PETSc also has other data structures and methods called index sets, and these methods act on PETSc vectors transform the data between the two views as required. These have also been implemented, in particular, a module was developed to transform vectors between the linearised- and three-vector data structures.

A code fragment showing the use of index sets to transform from the linearised vector, to the three vector is shown below.

SUBROUTINE xyzToComponents(linearVec,compVec,Nsize)
! converts a linearized vector of size 3*Nsize into an
! array of size 3, vectors of size N
! i.e. convert to a three-vector of x,y and z components
#include "finclude/petscsysdef.h"
#include "finclude/petscvecdef.h"
use petscVec
implicit none
Vec linearVec
Vec ,pointer :: compVec(:)
PetscInt Nsize
integer :: ierr, comm
IS xyz, comp
VecScatter scatter
comm=PETSC_COMM_WORLD
! x compVec first
CALL ISCreateStride(comm,Nsize,0,3,xyz,ierr)
CALL ISCreateStride(comm,Nsize,0,1,comp,ierr)
CALL VecScatterCreate(linearVec,xyz,compVec(1),comp, $
scatter,ierr)
CALL VecScatterBegin(scatter,linearVec,compVec(1), $
INSERT_VALUES,SCATTER_FORWARD,ierr)
CALL VecScatterEnd(scatter,linearVec,compVec(1), $
INSERT_VALUES,SCATTER_FORWARD,ierr)
CALL VecScatterDestroy(scatter,ierr)
CALL ISDestroy(xyz,ierr)
! don’t need to destroy comp, as it is the same
! for x,y and z

The fragment shows the translation for the X component, and this code is similar for the Y and Z components. Another subroutine has been written to perform the inverse translation. The RHS function performs the following tasks:

• Convert local linearised array Y to 3-d array m .
• Calculate the total energy htot for each element by calling function totalFieldVec .
• Calculate the derivative of the total magnetic field YDOT , which depend on htot and m .

The CVODE solver then uses YDOT to update Y . These three tasks have been implemented.

The MicroMag code makes use of advanced FORTRAN 90 features, such as pointers, which are supported by the PETSc library, to create the three-vectors, as shown in the code fragment below.

Vec term1,term2
Vec , pointer :: htot(:), mxyz(:), ydot_xyz(:)
CALL VecCreateMPI(comm,nmax_l,NMAX,term1,mpierr)
ithree=3
CALL VecDuplicateVecsF90(term1,ithree,htot,mpierr)

The function totalFieldVec performs the following tasks to calculate the :

• total force totfx for the de-magnetisation by calling function calcdemagVec .
• anisotropic contribution to the energy, hanis .
• demagnetisation contribution to the energy, hdemag .
• exchange contribution to the energy, hexch .
• external contribution to the energy, hext .

The total energy htot is then determined from the sum of the components for hanis , hdemag , hexch and hext . These tasks have all been implemented.