1 // Copyright (c) 2017-2026, Lawrence Livermore National Security, LLC and other CEED contributors.
2 // All Rights Reserved. See the top-level LICENSE and NOTICE files for details.
4 // SPDX-License-Identifier: BSD-2-Clause
10 …imple usage of libCEED to compute the surface area of a 3D body using matrix-free application of a…
13 // The example has no dependencies, and is designed to be self-contained.
16 // All libCEED objects use a Ceed device object constructed based on a command line argument (-ceed…
20 //     make ex2-surface [CEED_DIR=</path/to/libceed>]
22 // Sample runs:
24 //     ./ex2-surface
25 //     ./ex2-surface -ceed /cpu/self
26 //     ./ex2-surface -ceed /gpu/cuda
28 // Test in 1D-3D
29 //TESTARGS(name="1D User QFunction") -ceed {ceed_resource} -d 1 -t
30 //TESTARGS(name="2D User QFunction") -ceed {ceed_resource} -d 2 -t
31 //TESTARGS(name="3D User QFunction") -ceed {ceed_resource} -d 3 -t
32 //TESTARGS(name="1D Gallery QFunction") -ceed {ceed_resource} -d 1 -t -g
33 //TESTARGS(name="2D Gallery QFunction") -ceed {ceed_resource} -d 2 -t -g
34 //TESTARGS(name="3D Gallery QFunction") -ceed {ceed_resource} -d 3 -t -g
39 #include "ex2-surface.h"
62   CeedInt     prob_size   = -1;              // approximate problem size  in main()
69     if (!strcmp(argv[ia], "-h")) {  in main()
71     } else if (!strcmp(argv[ia], "-c") || !strcmp(argv[ia], "-ceed")) {  in main()
73     } else if (!strcmp(argv[ia], "-d")) {  in main()
76     } else if (!strcmp(argv[ia], "-m")) {  in main()
78     } else if (!strcmp(argv[ia], "-p")) {  in main()
80     } else if (!strcmp(argv[ia], "-q")) {  in main()
82     } else if (!strcmp(argv[ia], "-s")) {  in main()
84     } else if (!strcmp(argv[ia], "-b")) {  in main()
86     } else if (!strcmp(argv[ia], "-t")) {  in main()
88     } else if (!strcmp(argv[ia], "-g")) {  in main()
107     printf("  Ceed specification     [-c] : %s\n", ceed_spec);  in main()
108     printf("  Mesh dimension         [-d] : %" CeedInt_FMT "\n", dim);  in main()
109     printf("  Mesh degree            [-m] : %" CeedInt_FMT "\n", mesh_degree);  in main()
110     printf("  Solution degree        [-p] : %" CeedInt_FMT "\n", sol_degree);  in main()
111     printf("  Num. 1D quadrature pts [-q] : %" CeedInt_FMT "\n", num_qpts);  in main()
112     printf("  Approx. # unknowns     [-s] : %" CeedInt_FMT "\n", prob_size);  in main()
113     printf("  QFunction source       [-g] : %s\n", gallery ? "gallery" : "header");  in main()
115       printf("Test/quiet mode is %s\n", (test ? "ON" : "OFF (use -t to enable)"));  in main()
122   // Select appropriate backend and logical device based on the (-ceed) command line argument.  in main()
133   // Determine the mesh size based on the given approximate problem size.  in main()
237   // Create auxiliary solution-size vectors.  in main()
261   // Benchmark runs  in main()
264     printf(" Executing %d benchmarking runs...\n", benchmark);  in main()
287     printf("Surface area error         : % .14g\n", surface_area - exact_surface_area);  in main()
290     CeedScalar tol = (dim == 1 ? 10000. * CEED_EPSILON : dim == 2 ? 1E-1 : 1E-1);  in main()
292 …if (fabs(surface_area - exact_surface_area) > tol) printf("Surface area error         : % .14g\n",…  in main()
331       r--;  in GetCartesianMeshSize()
347     nd[d] = num_xyz[d] * (p - 1) + 1;  in BuildCartesianRestriction()
351   // elem:         0             1                 n-1  in BuildCartesianRestriction()
352   //           |---*-...-*---|---*-...-*---|- ... -|--...--|  in BuildCartesianRestriction()
353   // num_nodes:   0   1    p-1  p  p+1       2*p             n*p  in BuildCartesianRestriction()
369         g_nodes += (e_xyz[d] * (p - 1) + r_nodes % p) * g_nodes_stride;  in BuildCartesianRestriction()
394     nd[d] = num_xyz[d] * (p - 1) + 1;  in SetCartesianMeshCoords()
403   CeedLobattoQuadrature(p, nodes, NULL);  // nodes are in [-1,1]  in SetCartesianMeshCoords()
411       coords[gs_nodes + scalar_size * d] = ((d1d / (p - 1)) + nodes[d1d % (p - 1)]) / num_xyz[d];  in SetCartesianMeshCoords()
431     coords[i] = 0.5 + 1. / sqrt(3.) * sin((2. / 3.) * M_PI * (coords[i] - 0.5));  in TransformMeshCoords()