1. Lisandro Dalcin
  2. PetIGA

Commits

Lisandro Dalcin  committed 055cf57

Change calling convention of TSUserIFunction and TSUserIJacobian

  • Participants
  • Parent commits 5b18640
  • Branches default

Comments (0)

Files changed (5)

File demo/CahnHilliard2D.c

View file
 #undef  __FUNCT__
 #define __FUNCT__ "Residual"
 PetscErrorCode Residual(IGAPoint p,PetscReal dt,PetscReal shift,
-                        PetscReal t,const PetscScalar *U,const PetscScalar *V,
+                        PetscReal t,const PetscScalar *V,const PetscScalar *U,
                         PetscScalar *R,void *ctx)
 {
   AppCtx *user = (AppCtx *)ctx;
 
-  PetscScalar c,c_t;
+  PetscScalar c_t,c;
+  IGAPointInterpolate(p,0,V,&c_t);
   IGAPointInterpolate(p,0,U,&c);
-  IGAPointInterpolate(p,0,V,&c_t);
 
   PetscReal M,dM;
   Mobility(user,c,&M,&dM,NULL);
 #undef  __FUNCT__
 #define __FUNCT__ "Tangent"
 PetscErrorCode Tangent(IGAPoint p,PetscReal dt,PetscReal shift,
-                       PetscReal t,const PetscScalar *U,const PetscScalar *V,
+                       PetscReal t,const PetscScalar *V,const PetscScalar *U,
                        PetscScalar *K,void *ctx)
 {
   AppCtx *user = (AppCtx *)ctx;
 
-  PetscScalar c,c_t;
+  PetscScalar c_t,c;
+  IGAPointInterpolate(p,0,V,&c_t);
   IGAPointInterpolate(p,0,U,&c);
-  IGAPointInterpolate(p,0,V,&c_t);
 
   PetscReal M,dM,d2M;
   Mobility(user,c,&M,&dM,&d2M);

File demo/NavierStokesVMS.c

View file
 #undef  __FUNCT__
 #define __FUNCT__ "Residual"
 PetscErrorCode Residual(IGAPoint pnt,PetscReal dt,PetscReal shift,
-                        PetscReal t,const PetscScalar *U,const PetscScalar *V,
+                        PetscReal t,const PetscScalar *V,const PetscScalar *U,
                         PetscScalar *Re,void *ctx)
 {
   AppCtx *user = (AppCtx *)ctx;
   PetscReal nu = user->nu;
 
-  PetscScalar u[4],u_t[4];
+  PetscScalar u_t[4],u[4];
   PetscScalar grad_u[4][3];
   PetscScalar der2_u[4][3][3];
   IGAPointInterpolate(pnt,0,V,&u_t[0]);
 #undef  __FUNCT__
 #define __FUNCT__ "Tangent"
 PetscErrorCode Tangent(IGAPoint pnt,PetscReal dt,PetscReal shift,
-                       PetscReal t,const PetscScalar *U,const PetscScalar *V,
+                       PetscReal t,const PetscScalar *V,const PetscScalar *U,
                        PetscScalar *Ke,void *ctx)
 {
   AppCtx *user = (AppCtx *)ctx;

File demo/PatternFormation.c

View file
 #undef  __FUNCT__
 #define __FUNCT__ "Function"
 PetscErrorCode Function(IGAPoint point,PetscReal dt,PetscReal shift,
-                        PetscReal t,const PetscScalar *U,const PetscScalar *V,
+                        PetscReal t,const PetscScalar *V,const PetscScalar *U,
                         PetscScalar *Re,void *ctx)
 {
   AppCtx *user = (AppCtx *)ctx;
 #undef  __FUNCT__
 #define __FUNCT__ "Jacobian"
 PetscErrorCode Jacobian(IGAPoint point,PetscReal dt,PetscReal shift,
-                        PetscReal t,const PetscScalar *U,const PetscScalar *V,
+                        PetscReal t,const PetscScalar *V,const PetscScalar *U,
                         PetscScalar *Ke,void *ctx)
 {
   AppCtx *user = (AppCtx *)ctx;

File include/petiga.h

View file
 typedef PetscErrorCode (*IGAUserFunction)  (IGAPoint point,const PetscScalar *U,PetscScalar *F,void *ctx);
 typedef PetscErrorCode (*IGAUserJacobian)  (IGAPoint point,const PetscScalar *U,PetscScalar *J,void *ctx);
 typedef PetscErrorCode (*IGAUserIFunction) (IGAPoint point,PetscReal dt,PetscReal a,
-                                            PetscReal t,const PetscScalar *U,const PetscScalar *V,PetscScalar *F,
+                                            PetscReal t,const PetscScalar *V,const PetscScalar *U,PetscScalar *F,
                                             void *ctx);
 typedef PetscErrorCode (*IGAUserIJacobian) (IGAPoint point,PetscReal dt,PetscReal a,
-                                            PetscReal t,const PetscScalar *U,const PetscScalar *V,PetscScalar *J,
+                                            PetscReal t,const PetscScalar *V,const PetscScalar *U,PetscScalar *J,
                                             void *ctx);
 
 typedef struct _IGAUserOps *IGAUserOps;
 
 extern PetscErrorCode IGACreateTS(IGA iga,TS *ts);
 extern PetscErrorCode IGAFormIFunction(IGA iga,PetscReal dt,PetscReal a,
-                                       PetscReal t,Vec U,Vec V,Vec F,
+                                       PetscReal t,Vec V,Vec U,Vec F,
                                        IGAUserIFunction IFunction,void *ctx);
 extern PetscErrorCode IGAFormIJacobian(IGA iga,PetscReal dt,PetscReal a,
-                                       PetscReal t,Vec U,Vec V,Mat J,
+                                       PetscReal t,Vec V,Vec U,Mat J,
                                        IGAUserIJacobian IJacobian,void *ctx);
 
 /* ---------------------------------------------------------------- */

File src/petigats.c

View file
 #undef  __FUNCT__
 #define __FUNCT__ "IGAFormIFunction"
 PetscErrorCode IGAFormIFunction(IGA iga,PetscReal dt,PetscReal shift,
-                                PetscReal t,Vec vecU,Vec vecV,Vec vecF,
+                                PetscReal t,Vec vecV,Vec vecU,Vec vecF,
                                 IGAUserIFunction IFunction, void *ctx)
 {
+  Vec               localV;
   Vec               localU;
-  Vec               localV;
+  const PetscScalar *arrayV;
   const PetscScalar *arrayU;
-  const PetscScalar *arrayV;
   IGAElement        element;
   IGAPoint          point;
   PetscErrorCode    ierr;
   PetscFunctionBegin;
   PetscValidHeaderSpecific(iga,IGA_CLASSID,1);
-  PetscValidHeaderSpecific(vecU,VEC_CLASSID,2);
-  PetscValidHeaderSpecific(vecV,VEC_CLASSID,3);
-  PetscValidHeaderSpecific(vecF,VEC_CLASSID,4);
+  PetscValidHeaderSpecific(vecV,VEC_CLASSID,5);
+  PetscValidHeaderSpecific(vecU,VEC_CLASSID,6);
+  PetscValidHeaderSpecific(vecF,VEC_CLASSID,7);
 
   /* Clear global vector F */
   ierr = VecZeroEntries(vecF);CHKERRQ(ierr);
 
   /* Get local vectors U and V */
+  ierr = IGAGetLocalVec(iga,&localV);CHKERRQ(ierr);
   ierr = IGAGetLocalVec(iga,&localU);CHKERRQ(ierr);
-  ierr = IGAGetLocalVec(iga,&localV);CHKERRQ(ierr);
   /* Communicate global to local */
+  ierr = IGAGlobalToLocal(iga,vecV,localV);CHKERRQ(ierr);
   ierr = IGAGlobalToLocal(iga,vecU,localU);CHKERRQ(ierr);
-  ierr = IGAGlobalToLocal(iga,vecV,localV);CHKERRQ(ierr);
   /* Get array from the local vectors */
+  ierr = VecGetArrayRead(localV,&arrayV);CHKERRQ(ierr);
   ierr = VecGetArrayRead(localU,&arrayU);CHKERRQ(ierr);
-  ierr = VecGetArrayRead(localV,&arrayV);CHKERRQ(ierr);
 
   /* Element loop */
   ierr = IGAGetElement(iga,&element);CHKERRQ(ierr);
   ierr = IGAElementBegin(element);CHKERRQ(ierr);
   while (IGAElementNext(element)) {
-    PetscScalar *U, *V, *F;
+    PetscScalar *V, *U, *F;
+    ierr = IGAElementGetWorkVec(element,&V);CHKERRQ(ierr);
     ierr = IGAElementGetWorkVec(element,&U);CHKERRQ(ierr);
-    ierr = IGAElementGetWorkVec(element,&V);CHKERRQ(ierr);
     ierr = IGAElementGetWorkVec(element,&F);CHKERRQ(ierr);
     /* */
+    ierr = IGAElementGetValues(element,arrayV,V);CHKERRQ(ierr);
     ierr = IGAElementGetValues(element,arrayU,U);CHKERRQ(ierr);
-    ierr = IGAElementGetValues(element,arrayV,V);CHKERRQ(ierr);
     ierr = IGAElementFixValues(element,U);CHKERRQ(ierr);
     /* Quadrature loop */
     ierr = IGAElementGetPoint(element,&point);CHKERRQ(ierr);
     while (IGAPointNext(point)) {
       PetscScalar *R;
       ierr = IGAPointGetWorkVec(point,&R);CHKERRQ(ierr);
-      ierr = IFunction(point,dt,shift,t,U,V,R,ctx);CHKERRQ(ierr);
+      ierr = IFunction(point,dt,shift,t,V,U,R,ctx);CHKERRQ(ierr);
       ierr = IGAPointAddVec(point,R,F);CHKERRQ(ierr);
     }
     /* */
 #undef  __FUNCT__
 #define __FUNCT__ "IGAFormIJacobian"
 PetscErrorCode IGAFormIJacobian(IGA iga,PetscReal dt,PetscReal shift,
-                                PetscReal t,Vec vecU,Vec vecV,Mat matJ,
+                                PetscReal t,Vec vecV,Vec vecU,Mat matJ,
                                 IGAUserIJacobian IJacobian,void *ctx)
 {
+  Vec               localV;
   Vec               localU;
-  Vec               localV;
+  const PetscScalar *arrayV;
   const PetscScalar *arrayU;
-  const PetscScalar *arrayV;
   IGAElement        element;
   IGAPoint          point;
   PetscErrorCode    ierr;
   PetscFunctionBegin;
   PetscValidHeaderSpecific(iga,IGA_CLASSID,1);
-  PetscValidHeaderSpecific(vecU,VEC_CLASSID,3);
-  PetscValidHeaderSpecific(vecV,VEC_CLASSID,4);
-  PetscValidHeaderSpecific(matJ,MAT_CLASSID,6);
+  PetscValidHeaderSpecific(vecV,VEC_CLASSID,5);
+  PetscValidHeaderSpecific(vecU,VEC_CLASSID,6);
+  PetscValidHeaderSpecific(matJ,MAT_CLASSID,7);
 
   /* Clear global matrix J*/
   ierr = MatZeroEntries(matJ);CHKERRQ(ierr);
 
   /* Get local vectors U and V */
+  ierr = IGAGetLocalVec(iga,&localV);CHKERRQ(ierr);
   ierr = IGAGetLocalVec(iga,&localU);CHKERRQ(ierr);
-  ierr = IGAGetLocalVec(iga,&localV);CHKERRQ(ierr);
   /* Communicate global to local */
+  ierr = IGAGlobalToLocal(iga,vecV,localV);CHKERRQ(ierr);
   ierr = IGAGlobalToLocal(iga,vecU,localU);CHKERRQ(ierr);
-  ierr = IGAGlobalToLocal(iga,vecV,localV);CHKERRQ(ierr);
   /* Get array from the local vectors */
+  ierr = VecGetArrayRead(localV,&arrayV);CHKERRQ(ierr);
   ierr = VecGetArrayRead(localU,&arrayU);CHKERRQ(ierr);
-  ierr = VecGetArrayRead(localV,&arrayV);CHKERRQ(ierr);
 
   /* Element Loop */
   ierr = IGAGetElement(iga,&element);CHKERRQ(ierr);
   ierr = IGAElementBegin(element);CHKERRQ(ierr);
   while (IGAElementNext(element)) {
-    PetscScalar *U, *V, *J;
+    PetscScalar *V, *U, *J;
     ierr = IGAElementGetWorkVec(element,&U);CHKERRQ(ierr);
     ierr = IGAElementGetWorkVec(element,&V);CHKERRQ(ierr);
     ierr = IGAElementGetWorkMat(element,&J);CHKERRQ(ierr);
     /* */
+    ierr = IGAElementGetValues(element,arrayV,V);CHKERRQ(ierr);
     ierr = IGAElementGetValues(element,arrayU,U);CHKERRQ(ierr);
-    ierr = IGAElementGetValues(element,arrayV,V);CHKERRQ(ierr);
     ierr = IGAElementFixValues(element,U);CHKERRQ(ierr);
     /* Quadrature loop */
     ierr = IGAElementGetPoint(element,&point);CHKERRQ(ierr);
     while (IGAPointNext(point)) {
       PetscScalar *K;
       ierr = IGAPointGetWorkMat(point,&K);CHKERRQ(ierr);
-      ierr = IJacobian(point,dt,shift,t,U,V,K,ctx);CHKERRQ(ierr);
+      ierr = IJacobian(point,dt,shift,t,V,U,K,ctx);CHKERRQ(ierr);
       ierr = IGAPointAddMat(point,K,J);CHKERRQ(ierr);
     }
     /* */
   ierr = IGAElementEnd(element);CHKERRQ(ierr);
 
   /* Restore array to the vectors */
+  ierr = VecRestoreArrayRead(localV,&arrayV);CHKERRQ(ierr);
   ierr = VecRestoreArrayRead(localU,&arrayU);CHKERRQ(ierr);
-  ierr = VecRestoreArrayRead(localV,&arrayV);CHKERRQ(ierr);
   /* Restore local vectors U and V */
   ierr = IGARestoreLocalVec(iga,&localV);CHKERRQ(ierr);
   ierr = IGARestoreLocalVec(iga,&localU);CHKERRQ(ierr);
   PetscValidHeaderSpecific(F,VEC_CLASSID,5);
   PetscValidHeaderSpecific(iga,IGA_CLASSID,6);
   ierr = TSGetTimeStep(ts,&dt);CHKERRQ(ierr);
-  ierr = IGAFormIFunction(iga,dt,shift,t,U,V,F,
+  ierr = IGAFormIFunction(iga,dt,shift,t,V,U,F,
                           iga->userops->IFunction,
                           iga->userops->IFunCtx);CHKERRQ(ierr);
   PetscFunctionReturn(0);
   PetscValidPointer(m,8);
   PetscValidHeaderSpecific(iga,IGA_CLASSID,9);
   ierr = TSGetTimeStep(ts,&dt);CHKERRQ(ierr);
-  ierr = IGAFormIJacobian(iga,dt,shift,t,U,V,*P,
+  ierr = IGAFormIJacobian(iga,dt,shift,t,V,U,*P,
                           iga->userops->IJacobian,
                           iga->userops->IJacCtx);CHKERRQ(ierr);
   if (*J != * P) {