Skip to content

Commit

Permalink
suncomplextype
Browse files Browse the repository at this point in the history
  • Loading branch information
maggul committed Sep 28, 2024
1 parent 418f60b commit 9eb4bff
Show file tree
Hide file tree
Showing 23 changed files with 310 additions and 271 deletions.
14 changes: 7 additions & 7 deletions examples/arkode/CXX_parallel/ark_heat2D_lsrk_p.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -196,8 +196,8 @@ struct UserData
static int f(sunrealtype t, N_Vector u, N_Vector f, void* user_data);

// Spectral radius estimation routine
static int eig(sunrealtype t, N_Vector y, N_Vector fn, sunrealtype* lambdaR,
sunrealtype* lambdaI, void* user_data, N_Vector temp1,
static int eig(sunrealtype t, N_Vector y, N_Vector fn,
suncomplextype* lambda, void* user_data, N_Vector temp1,
N_Vector temp2, N_Vector temp3);

// -----------------------------------------------------------------------------
Expand Down Expand Up @@ -1121,17 +1121,17 @@ static int WaitRecv(UserData* udata)
}

// Spectral radius estimation routine
static int eig(sunrealtype t, N_Vector y, N_Vector fn, sunrealtype* lambdaR,
sunrealtype* lambdaI, void* user_data, N_Vector temp1,
static int eig(sunrealtype t, N_Vector y, N_Vector fn,
suncomplextype* lambda, void* user_data, N_Vector temp1,
N_Vector temp2, N_Vector temp3)
{
// Access problem data
UserData* udata = (UserData*)user_data;

// Fill in spectral radius value
*lambdaR = -SUN_RCONST(8.0) * max(udata->kx / udata->dx / udata->dx,
udata->ky / udata->dy / udata->dy);
*lambdaI = SUN_RCONST(0.0);
*lambda = - SUN_RCONST(8.0) * max(udata->kx / udata->dx / udata->dx,
udata->ky / udata->dy / udata->dy)
+ SUN_RCONST(0.0) * SUN_I;

// return with success
return 0;
Expand Down
14 changes: 7 additions & 7 deletions examples/arkode/CXX_serial/ark_heat2D_lsrk.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -140,8 +140,8 @@ struct UserData
static int f(sunrealtype t, N_Vector u, N_Vector f, void* user_data);

// Spectral radius estimation routine
static int eig(sunrealtype t, N_Vector y, N_Vector fn, sunrealtype* lambdaR,
sunrealtype* lambdaI, void* user_data, N_Vector temp1,
static int eig(sunrealtype t, N_Vector y, N_Vector fn,
suncomplextype* lambda, void* user_data, N_Vector temp1,
N_Vector temp2, N_Vector temp3);

// -----------------------------------------------------------------------------
Expand Down Expand Up @@ -498,17 +498,17 @@ static int f(sunrealtype t, N_Vector u, N_Vector f, void* user_data)
}

// Spectral radius estimation routine
static int eig(sunrealtype t, N_Vector y, N_Vector fn, sunrealtype* lambdaR,
sunrealtype* lambdaI, void* user_data, N_Vector temp1,
static int eig(sunrealtype t, N_Vector y, N_Vector fn,
suncomplextype* lambda, void* user_data, N_Vector temp1,
N_Vector temp2, N_Vector temp3)
{
// Access problem data
UserData* udata = (UserData*)user_data;

// Fill in spectral radius value
*lambdaR = -SUN_RCONST(8.0) * max(udata->kx / udata->dx / udata->dx,
udata->ky / udata->dy / udata->dy);
*lambdaI = SUN_RCONST(0.0);
*lambda = - SUN_RCONST(8.0) * max(udata->kx / udata->dx / udata->dx,
udata->ky / udata->dy / udata->dy)
+ SUN_RCONST(0.0) * SUN_I;

// return with success
return 0;
Expand Down
29 changes: 15 additions & 14 deletions examples/arkode/C_serial/ark_analytic_lsrk.c
Original file line number Diff line number Diff line change
Expand Up @@ -15,11 +15,11 @@
*
* The following is a simple example problem with analytical
* solution,
* dy/dt = lambda*y + 1/(1+t^2) - lambda*atan(t)
* dy/dt = gamma*y + 1/(1+t^2) - gamma*atan(t)
* for t in the interval [0.0, 10.0], with initial condition: y=0.
*
* The stiffness of the problem is directly proportional to the
* value of "lambda". The value of lambda should be negative to
* value of "gamma". The value of gamma should be negative to
* result in a well-posed ODE; for values with magnitude larger
* than 100 the problem becomes quite stiff.
*
Expand Down Expand Up @@ -58,8 +58,8 @@
static int f(sunrealtype t, N_Vector y, N_Vector ydot, void* user_data);

/* User-supplied Dominated Eigenvalue Called by the Solver */
static int dom_eig(sunrealtype t, N_Vector y, N_Vector fn, sunrealtype* lambdaR,
sunrealtype* lambdaI, void* user_data, N_Vector temp1,
static int dom_eig(sunrealtype t, N_Vector y, N_Vector fn,
suncomplextype* lambda, void* user_data, N_Vector temp1,
N_Vector temp2, N_Vector temp3);

/* Private function to check function return values */
Expand All @@ -82,7 +82,7 @@ int main(void)
sunindextype NEQ = 1; /* number of dependent vars. */
sunrealtype reltol = SUN_RCONST(1.0e-8); /* tolerances */
sunrealtype abstol = SUN_RCONST(1.0e-8);
sunrealtype lambda = SUN_RCONST(-1000000.0); /* stiffness parameter */
sunrealtype gamma = SUN_RCONST(-1000000.0); /* stiffness parameter */

/* general problem variables */
int flag; /* reusable error-checking flag */
Expand All @@ -98,7 +98,7 @@ int main(void)

/* Initial diagnostics output */
printf("\nAnalytical ODE test problem:\n");
printf(" lambda = %" GSYM "\n", lambda);
printf(" gamma = %" GSYM "\n", gamma);
printf(" reltol = %.1" ESYM "\n", reltol);
printf(" abstol = %.1" ESYM "\n\n", abstol);

Expand All @@ -115,7 +115,7 @@ int main(void)

/* Set routines */
flag = ARKodeSetUserData(arkode_mem,
(void*)&lambda); /* Pass lambda to user functions */
(void*)&gamma); /* Pass gamma to user functions */
if (check_flag(&flag, "ARKodeSetUserData", 1)) { return 1; }

/* Specify tolerances */
Expand Down Expand Up @@ -210,25 +210,26 @@ int main(void)
static int f(sunrealtype t, N_Vector y, N_Vector ydot, void* user_data)
{
sunrealtype* rdata = (sunrealtype*)user_data; /* cast user_data to sunrealtype */
sunrealtype lambda = rdata[0]; /* set shortcut for stiffness parameter */
sunrealtype gamma = rdata[0]; /* set shortcut for stiffness parameter */
sunrealtype u = N_VGetArrayPointer(y)[0]; /* access current solution value */

/* fill in the RHS function: "N_VGetArrayPointer" accesses the 0th entry of ydot */
N_VGetArrayPointer(ydot)[0] =
lambda * u + SUN_RCONST(1.0) / (SUN_RCONST(1.0) + t * t) - lambda * ATAN(t);
gamma * u + SUN_RCONST(1.0) / (SUN_RCONST(1.0) + t * t) - gamma * ATAN(t);

return 0; /* return with success */
}

/* dom_eig routine to estimate the dominated eigenvalue */
static int dom_eig(sunrealtype t, N_Vector y, N_Vector fn, sunrealtype* lambdaR,
sunrealtype* lambdaI, void* user_data, N_Vector temp1,
static int dom_eig(sunrealtype t, N_Vector y, N_Vector fn,
suncomplextype* lambda, void* user_data, N_Vector temp1,
N_Vector temp2, N_Vector temp3)
{
sunrealtype* rdata = (sunrealtype*)user_data; /* cast user_data to sunrealtype */
sunrealtype lambda = rdata[0]; /* set shortcut for stiffness parameter */
*lambdaR = lambda;
*lambdaI = 0.0;
sunrealtype gamma = rdata[0]; /* set shortcut for stiffness parameter */

*lambda = gamma + SUN_RCONST(0.0) * SUN_I;

return 0; /* return with success */
}

Expand Down
4 changes: 2 additions & 2 deletions examples/arkode/C_serial/ark_analytic_lsrk.out
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@

Analytical ODE test problem:
lambda = -1e+06
gamma = -1e+06
reltol = 1.0e-08
abstol = 1.0e-08

Expand Down Expand Up @@ -39,4 +39,4 @@ Max. num. of stages allowed = 200
Max. spectral radius = 1010000.00
Min. spectral radius = 1010000.00

ACCURACY at the final time = 1.01918e-13
ACCURACY at the final time = 1.01918e-13
54 changes: 27 additions & 27 deletions examples/arkode/C_serial/ark_analytic_lsrk_varjac.c
Original file line number Diff line number Diff line change
Expand Up @@ -15,20 +15,20 @@
*
* The following is a simple example problem that solves the ODE
*
* dy/dt = (lambda - alpha*cos((10 - t)/10*pi)*y + 1/(1+t^2)
* - (lambda - alpha*cos((10 - t)/10*pi)*atan(t),
* dy/dt = (gamma - theta*cos((10 - t)/10*pi)*y + 1/(1+t^2)
* - (gamma - theta*cos((10 - t)/10*pi)*atan(t),
*
* for t in the interval [0.0, 10.0], with an initial condition: y=0.
* This initial value problem has the analytical solution
*
* y(t) = arctan(t).
*
* The stiffness of the problem depends on both lambda and alpha together.
* While lambda determines the center of the stiffness parameter,
* the value of alpha determines the radius of the interval in which
* The stiffness of the problem depends on both gamma and theta together.
* While gamma determines the center of the stiffness parameter,
* the value of theta determines the radius of the interval in which
* the stiffness parameter lies.
*
* The value of "lambda - alpha*cos((10 - t)/10*pi)" should
* The value of "gamma - theta*cos((10 - t)/10*pi)" should
* be negative to result in a well-posed ODE; for values with magnitude
* larger than 100 the problem becomes quite stiff.
*
Expand Down Expand Up @@ -73,8 +73,8 @@
static int f(sunrealtype t, N_Vector y, N_Vector ydot, void* user_data);

/* User-supplied Dominated Eigenvalue Called by the Solver */
static int dom_eig(sunrealtype t, N_Vector y, N_Vector fn, sunrealtype* lambdaR,
sunrealtype* lambdaI, void* user_data, N_Vector temp1,
static int dom_eig(sunrealtype t, N_Vector y, N_Vector fn,
suncomplextype* lambda, void* user_data, N_Vector temp1,
N_Vector temp2, N_Vector temp3);

/* Private function to check function return values */
Expand All @@ -97,11 +97,11 @@ int main(void)
sunindextype NEQ = 1; /* number of dependent vars. */
sunrealtype reltol = SUN_RCONST(1.0e-8); /* tolerances */
sunrealtype abstol = SUN_RCONST(1.0e-8);
sunrealtype lambda = SUN_RCONST(-1.0e+6); /* stiffness parameter 1*/
sunrealtype alpha = SUN_RCONST(1.0e+2); /* stiffness parameter 2*/
sunrealtype gamma = SUN_RCONST(-1.0e+6); /* stiffness parameter 1*/
sunrealtype theta = SUN_RCONST(1.0e+2); /* stiffness parameter 2*/
sunrealtype UserData[2];
UserData[0] = lambda;
UserData[1] = alpha;
UserData[0] = gamma;
UserData[1] = theta;

/* general problem variables */
int flag; /* reusable error-checking flag */
Expand All @@ -118,9 +118,9 @@ int main(void)
/* Initial diagnostics output */
printf("\nAnalytical ODE test problem with a variable Jacobian:");
printf("\nThe stiffness of the problem is directly proportional to");
printf("\n\"lambda - alpha*cos((10 - t)/10*pi)\"\n\n");
printf(" lambda = %" GSYM "\n", lambda);
printf(" alpha = %" GSYM "\n", alpha);
printf("\n\"gamma - theta*cos((10 - t)/10*pi)\"\n\n");
printf(" gamma = %" GSYM "\n", gamma);
printf(" theta = %" GSYM "\n", theta);
printf(" reltol = %.1" ESYM "\n", reltol);
printf(" abstol = %.1" ESYM "\n\n", abstol);

Expand All @@ -137,7 +137,7 @@ int main(void)

/* Set routines */
flag = ARKodeSetUserData(arkode_mem,
(void*)&UserData); /* Pass lambda to user functions */
(void*)&UserData); /* Pass gamma to user functions */
if (check_flag(&flag, "ARKodeSetUserData", 1)) { return 1; }

/* Specify tolerances */
Expand Down Expand Up @@ -231,31 +231,31 @@ int main(void)
static int f(sunrealtype t, N_Vector y, N_Vector ydot, void* user_data)
{
sunrealtype* rdata = (sunrealtype*)user_data; /* cast user_data to sunrealtype */
sunrealtype lambda = rdata[0]; /* set shortcut for stiffness parameter 1 */
sunrealtype alpha = rdata[1]; /* set shortcut for stiffness parameter 2 */
sunrealtype gamma = rdata[0]; /* set shortcut for stiffness parameter 1 */
sunrealtype theta = rdata[1]; /* set shortcut for stiffness parameter 2 */
sunrealtype u = N_VGetArrayPointer(y)[0]; /* access current solution value */

/* fill in the RHS function: "N_VGetArrayPointer" accesses the 0th entry of ydot */
N_VGetArrayPointer(ydot)[0] =
(lambda - alpha * COS((10 - t) / 10 * ACOS(-1))) * u +
(gamma - theta * COS((10 - t) / 10 * ACOS(-1))) * u +
SUN_RCONST(1.0) / (SUN_RCONST(1.0) + t * t) -
(lambda - alpha * COS((10 - t) / 10 * ACOS(-1))) * ATAN(t);
(gamma - theta * COS((10 - t) / 10 * ACOS(-1))) * ATAN(t);

return 0; /* return with success */
}

/* dom_eig routine to estimate the dominated eigenvalue */
static int dom_eig(sunrealtype t, N_Vector y, N_Vector fn, sunrealtype* lambdaR,
sunrealtype* lambdaI, void* user_data, N_Vector temp1,
static int dom_eig(sunrealtype t, N_Vector y, N_Vector fn,
suncomplextype* lambda, void* user_data, N_Vector temp1,
N_Vector temp2, N_Vector temp3)
{
sunrealtype* rdata = (sunrealtype*)user_data; /* cast user_data to sunrealtype */
sunrealtype lambda = rdata[0]; /* set shortcut for stiffness parameter 1 */
sunrealtype alpha = rdata[1]; /* set shortcut for stiffness parameter 2 */
*lambdaR = (lambda - alpha * COS((10 - t) / 10 *
ACOS(-1))); /* access current solution value */
*lambdaI = 0.0;
sunrealtype gamma = rdata[0]; /* set shortcut for stiffness parameter 1 */
sunrealtype theta = rdata[1]; /* set shortcut for stiffness parameter 2 */

*lambda = (gamma - theta * COS((10 - t) / 10 * ACOS(-1)))
+ SUN_RCONST(0.0) * SUN_I;

return 0; /* return with success */
}

Expand Down
6 changes: 3 additions & 3 deletions examples/arkode/C_serial/ark_analytic_lsrk_varjac.out
Original file line number Diff line number Diff line change
@@ -1,10 +1,10 @@

Analytical ODE test problem with a variable Jacobian:
The stiffness of the problem is directly proportional to
"lambda - alpha*cos((10 - t)/10*pi)"
"gamma - theta*cos((10 - t)/10*pi)"

lambda = -1e+06
alpha = 100
gamma = -1e+06
theta = 100
reltol = 1.0e-08
abstol = 1.0e-08

Expand Down
12 changes: 6 additions & 6 deletions examples/sunlinsol/dense/test_sunlinsol_dense.c
Original file line number Diff line number Diff line change
Expand Up @@ -45,7 +45,7 @@ int main(int argc, char* argv[])
int fails = 0; /* counter for test failures */
sunindextype cols, rows; /* matrix columns, rows */
SUNLinearSolver LS; /* solver object */
SUNMatrix A, B, I; /* test matrices */
SUNMatrix A, B, K; /* test matrices */
N_Vector x, y, b; /* test vectors */
int print_timing;
int print_on_fail;
Expand Down Expand Up @@ -86,7 +86,7 @@ int main(int argc, char* argv[])
/* Create matrices and vectors */
A = SUNDenseMatrix(rows, cols, sunctx);
B = SUNDenseMatrix(rows, cols, sunctx);
I = SUNDenseMatrix(rows, cols, sunctx);
K = SUNDenseMatrix(rows, cols, sunctx);
x = N_VNew_Serial(cols, sunctx);
y = N_VNew_Serial(cols, sunctx);
b = N_VNew_Serial(cols, sunctx);
Expand All @@ -105,7 +105,7 @@ int main(int argc, char* argv[])
j = cols - 1;
for (k = 0; k < rows; k++)
{
colj = SUNDenseMatrix_Column(I, j);
colj = SUNDenseMatrix_Column(K, j);
colj[k] = 1;
j = j - 1;
}
Expand All @@ -116,7 +116,7 @@ int main(int argc, char* argv[])
for (j = 0; j < cols; j++)
{
colj = SUNDenseMatrix_Column(A, j);
colIj = SUNDenseMatrix_Column(I, j);
colIj = SUNDenseMatrix_Column(K, j);
colj[k] = colj[k] + colIj[k];
}
}
Expand All @@ -141,7 +141,7 @@ int main(int argc, char* argv[])
/* Free matrices and vectors */
SUNMatDestroy(A);
SUNMatDestroy(B);
SUNMatDestroy(I);
SUNMatDestroy(K);
N_VDestroy(x);
N_VDestroy(y);
N_VDestroy(b);
Expand Down Expand Up @@ -184,7 +184,7 @@ int main(int argc, char* argv[])
SUNLinSolFree(LS);
SUNMatDestroy(A);
SUNMatDestroy(B);
SUNMatDestroy(I);
SUNMatDestroy(K);
N_VDestroy(x);
N_VDestroy(y);
N_VDestroy(b);
Expand Down
Loading

0 comments on commit 9eb4bff

Please sign in to comment.