From f9eaf6bb83256b84fac7fa8c2d0c0e303035f719 Mon Sep 17 00:00:00 2001 From: allen Date: Wed, 5 Jan 2000 14:06:33 +0000 Subject: Tidying and adding error codes. Also made each of the registered solvers return an integer rather than a void so that the error codes can be passed back Standard elliptic error codes are now in CactusElliptic/EllBase/src/EllBase.h and so far are ELL_SUCCESS ELL_NOCONVERGENCE ELL_NOSOLVER git-svn-id: http://svn.cactuscode.org/arrangements/CactusElliptic/EllSOR/trunk@45 fa3da13c-9f13-4301-a575-cf5b8c5e1907 --- src/sor_flat.F | 18 ++++++----- src/sor_wrapper.c | 92 +++++++++++++++++++++++++++++++++++++++---------------- 2 files changed, 77 insertions(+), 33 deletions(-) diff --git a/src/sor_flat.F b/src/sor_flat.F index e36514a..ac78d52 100644 --- a/src/sor_flat.F +++ b/src/sor_flat.F @@ -11,6 +11,7 @@ #include "cctk_arguments.h" #include "cctk_parameters.h" +#include "CactusElliptic/EllBase/src/EllBase.h" /*@@ @routine sor @@ -23,7 +24,8 @@ @enddesc @@*/ - subroutine sor_flat_core3d(_CCTK_FARGUMENTS, + subroutine sor_flat_core3d( + & ierr,_CCTK_FARGUMENTS, $ Mlinear_lsh,Mlinear, $ Nsource_lsh,Nsource, $ var,var_idx, @@ -33,8 +35,9 @@ _DECLARE_CCTK_FARGUMENTS DECLARE_CCTK_PARAMETERS - INTEGER CCTK_Equals - + DECLARE_CCTK_FUNCTIONS + + integer ierr INTEGER Mlinear_lsh(3) CCTK_REAL Mlinear(Mlinear_lsh(1),Mlinear_lsh(2),Mlinear_lsh(3)) @@ -87,11 +90,13 @@ c Coeeficients for the solver: 19 point stencil... logical cheb, const, none, verb integer Mlinear_storage,Nsource_storage - INTEGER sum_handle,ierr + INTEGER sum_handle c stencil size INTEGER sw(3) + ierr = ELL_SUCCESS + tol = AbsTol(1) c Get the reduction handel for the sum operation @@ -248,9 +253,8 @@ c Apply octant Symmetries enddo - write (*,*) "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" - write (*,*) "!! WARNING: SOR SOLVER DID NOT CONVERGE !!" - write (*,*) "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" + call CCTK_INFO("SOR solver did not converge") + ierr = ELL_NOCONVERGENCE 123 continue diff --git a/src/sor_wrapper.c b/src/sor_wrapper.c index a0b3eac..d1ce7fa 100644 --- a/src/sor_wrapper.c +++ b/src/sor_wrapper.c @@ -22,24 +22,47 @@ #include "cctk_parameters.h" #include "cctk_FortranString.h" -void FORTRAN_NAME(sor_confmetric_core3d)(_CCTK_C2F_PROTO(GH), - int *, CCTK_REAL *, - int *, CCTK_REAL *, - CCTK_REAL *, - CCTK_REAL *,CCTK_REAL *,CCTK_REAL *, - CCTK_REAL *,CCTK_REAL *,CCTK_REAL *, - CCTK_REAL *, int *, CCTK_REAL *, CCTK_REAL *); - -void FORTRAN_NAME(sor_flat_core3d)(_CCTK_C2F_PROTO(GH), - int *, CCTK_REAL *, - int *, CCTK_REAL *, - CCTK_REAL *, int *, CCTK_REAL *, CCTK_REAL *); +void FORTRAN_NAME(sor_confmetric_core3d) + (_CCTK_C2F_PROTO(GH), + int *, + CCTK_REAL *, + int *, + CCTK_REAL *, + CCTK_REAL *, + CCTK_REAL *, + CCTK_REAL *, + CCTK_REAL *, + CCTK_REAL *, + CCTK_REAL *, + CCTK_REAL *, + CCTK_REAL *, + int *, + CCTK_REAL *, + CCTK_REAL *); + +void FORTRAN_NAME(sor_flat_core3d) + (int *ierr, + _CCTK_C2F_PROTO(GH), + int *, + CCTK_REAL *, + int *, + CCTK_REAL *, + CCTK_REAL *, + int *, + CCTK_REAL *, + CCTK_REAL *); /* We pass in the arguments that are neccessary for this class of elliptic eq. this solver is intended to solve. See ./CactusElliptic/EllBase/src/ for the classes of elliptic eq. */ -void sor_confmetric(cGH *GH, int *MetricPsiI, int FieldIndex, - int MIndex, int NIndex, CCTK_REAL *AbsTol,CCTK_REAL *RelTol) { +void sor_confmetric(cGH *GH, + int *MetricPsiI, + int FieldIndex, + int MIndex, + int NIndex, + CCTK_REAL *AbsTol, + CCTK_REAL *RelTol) +{ CCTK_REAL *gxx=NULL, *gxy=NULL, *gxz=NULL; CCTK_REAL *gyy=NULL, *gyz=NULL, *gzz=NULL; @@ -58,6 +81,7 @@ void sor_confmetric(cGH *GH, int *MetricPsiI, int FieldIndex, /* derive the metric data pointer from the index array. Note the ordering. Also get datapointers to the field to solve for. All of these are mandatory */ + gxx = (CCTK_REAL*) CCTK_VarDataPtrI(GH, 0, MetricPsiI[0]); gxy = (CCTK_REAL*) CCTK_VarDataPtrI(GH, 0, MetricPsiI[1]); gxz = (CCTK_REAL*) CCTK_VarDataPtrI(GH, 0, MetricPsiI[2]); @@ -69,8 +93,10 @@ void sor_confmetric(cGH *GH, int *MetricPsiI, int FieldIndex, Field = (CCTK_REAL*) CCTK_VarDataPtrI(GH,0,FieldIndex); if ((!gxx)||(!gxy)||(!gxz)||(!gyy)||(!gyz)||(!gzz)||(!psi)||(!Field)) + { CCTK_WARN(0,"SOR_WRAPPER: One of the metric data fields, or the GF to solve could not be found!"); - + } + /* derive the data pointer for the fields. the M/N fields are not allocated (better: are of size 1), if the passed index is negative, or we get back an empty GF of size 1 */ @@ -85,7 +111,8 @@ void sor_confmetric(cGH *GH, int *MetricPsiI, int FieldIndex, if (GH->cctk_dim>3) CCTK_WARN(0,"This elliptic solver implementation does not do dimension>3!"); - for (i=0;icctk_dim;i++) { + for (i=0;icctk_dim;i++) + { if((MIndex<0)) Mlinear_lsh[i]=1; else Mlinear_lsh[i]=GH->cctk_lsh[i]; if((NIndex<0)) Nsource_lsh[i]=1; @@ -103,20 +130,24 @@ void sor_confmetric(cGH *GH, int *MetricPsiI, int FieldIndex, } -void sor_flat(cGH *GH, int FieldIndex, int MIndex, - int NIndex, CCTK_REAL *AbsTol, CCTK_REAL *RelTol) { - +int sor_flat(cGH *GH, + int FieldIndex, + int MIndex, + int NIndex, + CCTK_REAL *AbsTol, + CCTK_REAL *RelTol) +{ + int ierr; + int retval=0; CCTK_REAL *Mlinear=NULL, *Nsources=NULL; CCTK_REAL *Field=NULL; CCTK_REAL tolerance; int i; - int toltype; - int Mlinear_lsh[3], Nsource_lsh[3]; - int retcode; Field = (CCTK_REAL*) CCTK_VarDataPtrI(GH,0,FieldIndex); + if (MIndex>0) Mlinear = (CCTK_REAL*) CCTK_VarDataPtrI(GH,0,MIndex); if (NIndex>0) Nsources = (CCTK_REAL*) CCTK_VarDataPtrI(GH,0,NIndex); @@ -131,10 +162,19 @@ void sor_flat(cGH *GH, int FieldIndex, int MIndex, } /* call the fortran routine */ - FORTRAN_NAME(sor_flat_core3d)(_PASS_CCTK_C2F(GH), - Mlinear_lsh, Mlinear, - Nsource_lsh, Nsources, - Field, &FieldIndex, AbsTol, RelTol); + FORTRAN_NAME(sor_flat_core3d) + (&ierr, + _PASS_CCTK_C2F(GH), + Mlinear_lsh, + Mlinear, + Nsource_lsh, + Nsources, + Field, + &FieldIndex, + AbsTol, + RelTol); + + return ierr; } -- cgit v1.2.3