xref: /libCEED/examples/solids/problems/neo-hookean.c (revision 5754ecac3b7d1ff97b39b25dc78c06350f2c900d)
1*5754ecacSJeremy L Thompson #include <ceed.h>
2*5754ecacSJeremy L Thompson #include <petsc.h>
3*5754ecacSJeremy L Thompson #include "../problems/neo-hookean.h"
4*5754ecacSJeremy L Thompson 
5*5754ecacSJeremy L Thompson // Build libCEED context object
6*5754ecacSJeremy L Thompson PetscErrorCode PhysicsContext_NH(MPI_Comm comm, Ceed ceed, Units *units,
7*5754ecacSJeremy L Thompson                                  CeedQFunctionContext *ctx) {
8*5754ecacSJeremy L Thompson   PetscErrorCode ierr;
9*5754ecacSJeremy L Thompson   Physics_NH phys;
10*5754ecacSJeremy L Thompson 
11*5754ecacSJeremy L Thompson   PetscFunctionBegin;
12*5754ecacSJeremy L Thompson 
13*5754ecacSJeremy L Thompson   ierr = PetscMalloc1(1, units); CHKERRQ(ierr);
14*5754ecacSJeremy L Thompson   ierr = PetscMalloc1(1, &phys); CHKERRQ(ierr);
15*5754ecacSJeremy L Thompson   ierr = ProcessPhysics_NH(comm, phys, *units); CHKERRQ(ierr);
16*5754ecacSJeremy L Thompson   CeedQFunctionContextCreate(ceed, ctx);
17*5754ecacSJeremy L Thompson   CeedQFunctionContextSetData(*ctx, CEED_MEM_HOST, CEED_COPY_VALUES,
18*5754ecacSJeremy L Thompson                               sizeof(*phys), phys);
19*5754ecacSJeremy L Thompson   ierr = PetscFree(phys); CHKERRQ(ierr);
20*5754ecacSJeremy L Thompson 
21*5754ecacSJeremy L Thompson   PetscFunctionReturn(0);
22*5754ecacSJeremy L Thompson }
23*5754ecacSJeremy L Thompson 
24*5754ecacSJeremy L Thompson // Build libCEED smoother context object
25*5754ecacSJeremy L Thompson PetscErrorCode PhysicsSmootherContext_NH(MPI_Comm comm, Ceed ceed,
26*5754ecacSJeremy L Thompson     CeedQFunctionContext ctx, CeedQFunctionContext *ctx_smoother) {
27*5754ecacSJeremy L Thompson   PetscErrorCode ierr;
28*5754ecacSJeremy L Thompson   PetscScalar nu_smoother = 0;
29*5754ecacSJeremy L Thompson   PetscBool nu_flag = PETSC_FALSE;
30*5754ecacSJeremy L Thompson   Physics_NH phys, phys_smoother;
31*5754ecacSJeremy L Thompson 
32*5754ecacSJeremy L Thompson   PetscFunctionBegin;
33*5754ecacSJeremy L Thompson 
34*5754ecacSJeremy L Thompson   ierr = PetscOptionsBegin(comm, NULL,
35*5754ecacSJeremy L Thompson                            "Neo-Hookean physical parameters for smoother", NULL);
36*5754ecacSJeremy L Thompson   CHKERRQ(ierr);
37*5754ecacSJeremy L Thompson 
38*5754ecacSJeremy L Thompson   ierr = PetscOptionsScalar("-nu_smoother", "Poisson's ratio for smoother",
39*5754ecacSJeremy L Thompson                             NULL, nu_smoother, &nu_smoother, &nu_flag);
40*5754ecacSJeremy L Thompson   CHKERRQ(ierr);
41*5754ecacSJeremy L Thompson 
42*5754ecacSJeremy L Thompson   ierr = PetscOptionsEnd(); CHKERRQ(ierr); // End of setting Physics
43*5754ecacSJeremy L Thompson 
44*5754ecacSJeremy L Thompson   if (nu_flag) {
45*5754ecacSJeremy L Thompson     // Copy context
46*5754ecacSJeremy L Thompson     CeedQFunctionContextGetData(ctx, CEED_MEM_HOST, &phys);
47*5754ecacSJeremy L Thompson     ierr = PetscMalloc1(1, &phys_smoother); CHKERRQ(ierr);
48*5754ecacSJeremy L Thompson     ierr = PetscMemcpy(phys_smoother, phys, sizeof(*phys)); CHKERRQ(ierr);
49*5754ecacSJeremy L Thompson     CeedQFunctionContextRestoreData(ctx, &phys);
50*5754ecacSJeremy L Thompson     // Create smoother context
51*5754ecacSJeremy L Thompson     CeedQFunctionContextCreate(ceed, ctx_smoother);
52*5754ecacSJeremy L Thompson     phys_smoother->nu = nu_smoother;
53*5754ecacSJeremy L Thompson     CeedQFunctionContextSetData(*ctx_smoother, CEED_MEM_HOST, CEED_COPY_VALUES,
54*5754ecacSJeremy L Thompson                                 sizeof(*phys_smoother), phys_smoother);
55*5754ecacSJeremy L Thompson     ierr = PetscFree(phys_smoother); CHKERRQ(ierr);
56*5754ecacSJeremy L Thompson   } else {
57*5754ecacSJeremy L Thompson     *ctx_smoother = NULL;
58*5754ecacSJeremy L Thompson   }
59*5754ecacSJeremy L Thompson 
60*5754ecacSJeremy L Thompson   PetscFunctionReturn(0);
61*5754ecacSJeremy L Thompson }
62*5754ecacSJeremy L Thompson 
63*5754ecacSJeremy L Thompson // Process physics options - Neo-Hookean
64*5754ecacSJeremy L Thompson PetscErrorCode ProcessPhysics_NH(MPI_Comm comm, Physics_NH phys, Units units) {
65*5754ecacSJeremy L Thompson   PetscErrorCode ierr;
66*5754ecacSJeremy L Thompson   PetscBool nu_flag = PETSC_FALSE;
67*5754ecacSJeremy L Thompson   PetscBool Young_flag = PETSC_FALSE;
68*5754ecacSJeremy L Thompson   phys->nu = 0;
69*5754ecacSJeremy L Thompson   phys->E = 0;
70*5754ecacSJeremy L Thompson   units->meter     = 1;        // 1 meter in scaled length units
71*5754ecacSJeremy L Thompson   units->second    = 1;        // 1 second in scaled time units
72*5754ecacSJeremy L Thompson   units->kilogram  = 1;        // 1 kilogram in scaled mass units
73*5754ecacSJeremy L Thompson 
74*5754ecacSJeremy L Thompson   PetscFunctionBeginUser;
75*5754ecacSJeremy L Thompson 
76*5754ecacSJeremy L Thompson   ierr = PetscOptionsBegin(comm, NULL, "Neo-Hookean physical parameters", NULL);
77*5754ecacSJeremy L Thompson   CHKERRQ(ierr);
78*5754ecacSJeremy L Thompson 
79*5754ecacSJeremy L Thompson   ierr = PetscOptionsScalar("-nu", "Poisson's ratio", NULL, phys->nu, &phys->nu,
80*5754ecacSJeremy L Thompson                             &nu_flag); CHKERRQ(ierr);
81*5754ecacSJeremy L Thompson 
82*5754ecacSJeremy L Thompson   ierr = PetscOptionsScalar("-E", "Young's Modulus", NULL, phys->E, &phys->E,
83*5754ecacSJeremy L Thompson                             &Young_flag); CHKERRQ(ierr);
84*5754ecacSJeremy L Thompson 
85*5754ecacSJeremy L Thompson   ierr = PetscOptionsScalar("-units_meter", "1 meter in scaled length units",
86*5754ecacSJeremy L Thompson                             NULL, units->meter, &units->meter, NULL);
87*5754ecacSJeremy L Thompson   CHKERRQ(ierr);
88*5754ecacSJeremy L Thompson   units->meter = fabs(units->meter);
89*5754ecacSJeremy L Thompson 
90*5754ecacSJeremy L Thompson   ierr = PetscOptionsScalar("-units_second", "1 second in scaled time units",
91*5754ecacSJeremy L Thompson                             NULL, units->second, &units->second, NULL);
92*5754ecacSJeremy L Thompson   CHKERRQ(ierr);
93*5754ecacSJeremy L Thompson   units->second = fabs(units->second);
94*5754ecacSJeremy L Thompson 
95*5754ecacSJeremy L Thompson   ierr = PetscOptionsScalar("-units_kilogram", "1 kilogram in scaled mass units",
96*5754ecacSJeremy L Thompson                             NULL, units->kilogram, &units->kilogram, NULL);
97*5754ecacSJeremy L Thompson   CHKERRQ(ierr);
98*5754ecacSJeremy L Thompson   units->kilogram = fabs(units->kilogram);
99*5754ecacSJeremy L Thompson 
100*5754ecacSJeremy L Thompson   ierr = PetscOptionsEnd(); CHKERRQ(ierr); // End of setting Physics
101*5754ecacSJeremy L Thompson 
102*5754ecacSJeremy L Thompson   // Check for all required options to be set
103*5754ecacSJeremy L Thompson   if (!nu_flag) {
104*5754ecacSJeremy L Thompson     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "-nu option needed");
105*5754ecacSJeremy L Thompson   }
106*5754ecacSJeremy L Thompson   if (!Young_flag) {
107*5754ecacSJeremy L Thompson     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "-E option needed");
108*5754ecacSJeremy L Thompson   }
109*5754ecacSJeremy L Thompson 
110*5754ecacSJeremy L Thompson   // Define derived units
111*5754ecacSJeremy L Thompson   units->Pascal = units->kilogram / (units->meter * PetscSqr(units->second));
112*5754ecacSJeremy L Thompson 
113*5754ecacSJeremy L Thompson   // Scale E to Pa
114*5754ecacSJeremy L Thompson   phys->E *= units->Pascal;
115*5754ecacSJeremy L Thompson 
116*5754ecacSJeremy L Thompson   PetscFunctionReturn(0);
117*5754ecacSJeremy L Thompson };