Compare commits
6
Commits
| Author | SHA1 | Date | |
|---|---|---|---|
 Tzanio
|
967908029f | ||
|
AMCBRIDGE\ahudym
|
8560248ea4 | ||
|
AMCBRIDGE\ahudym
|
fa85cf2471 | ||
|
AMCBRIDGE\ahudym
|
8b34881d46 | ||
|
AMCBRIDGE\ahudym
|
de0cad550d | ||
|
AMCBRIDGE\ahudym
|
1f7f84d7a2 |
@@ -175,6 +175,10 @@ else()
|
|||||||
set(MFEM_DEBUG OFF)
|
set(MFEM_DEBUG OFF)
|
||||||
endif()
|
endif()
|
||||||
|
|
||||||
|
if (WIN32)
|
||||||
|
add_definitions(-D_USE_MATH_DEFINES)
|
||||||
|
endif()
|
||||||
|
|
||||||
# MPI -> hypre; PETSc (optional)
|
# MPI -> hypre; PETSc (optional)
|
||||||
if (MFEM_USE_MPI)
|
if (MFEM_USE_MPI)
|
||||||
find_package(MPI REQUIRED)
|
find_package(MPI REQUIRED)
|
||||||
|
|||||||
+1
-1
@@ -92,7 +92,7 @@ class PMLDiagMatrixCoefficient : public VectorCoefficient
|
|||||||
{
|
{
|
||||||
private:
|
private:
|
||||||
CartesianPML * pml = nullptr;
|
CartesianPML * pml = nullptr;
|
||||||
void (*Function)(const Vector &, CartesianPML * , Vector &);
|
void (*Function)(const Vector &, CartesianPML *, Vector &);
|
||||||
public:
|
public:
|
||||||
PMLDiagMatrixCoefficient(int dim, void(*F)(const Vector &, CartesianPML *,
|
PMLDiagMatrixCoefficient(int dim, void(*F)(const Vector &, CartesianPML *,
|
||||||
Vector &),
|
Vector &),
|
||||||
|
|||||||
+1
-1
@@ -92,7 +92,7 @@ class PMLDiagMatrixCoefficient : public VectorCoefficient
|
|||||||
{
|
{
|
||||||
private:
|
private:
|
||||||
CartesianPML * pml = nullptr;
|
CartesianPML * pml = nullptr;
|
||||||
void (*Function)(const Vector &, CartesianPML * , Vector &);
|
void (*Function)(const Vector &, CartesianPML *, Vector &);
|
||||||
public:
|
public:
|
||||||
PMLDiagMatrixCoefficient(int dim, void(*F)(const Vector &, CartesianPML *,
|
PMLDiagMatrixCoefficient(int dim, void(*F)(const Vector &, CartesianPML *,
|
||||||
Vector &),
|
Vector &),
|
||||||
|
|||||||
@@ -17,14 +17,14 @@ namespace mfem
|
|||||||
{
|
{
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EAConvectionAssemble1D(const int NE,
|
void EAConvectionAssemble1D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata,
|
Vector &eadata,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -69,14 +69,14 @@ static void EAConvectionAssemble1D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EAConvectionAssemble2D(const int NE,
|
void EAConvectionAssemble2D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata,
|
Vector &eadata,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -146,14 +146,14 @@ static void EAConvectionAssemble2D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EAConvectionAssemble3D(const int NE,
|
void EAConvectionAssemble3D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata,
|
Vector &eadata,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
|
|||||||
@@ -21,13 +21,13 @@ namespace mfem
|
|||||||
// PA Convection Integrator
|
// PA Convection Integrator
|
||||||
|
|
||||||
// PA Convection Assemble 2D kernel
|
// PA Convection Assemble 2D kernel
|
||||||
static void PAConvectionSetup2D(const int NQ,
|
void PAConvectionSetup2D(const int NQ,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
const Vector &vel,
|
const Vector &vel,
|
||||||
const double alpha,
|
const double alpha,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
constexpr int DIM = 2;
|
constexpr int DIM = 2;
|
||||||
|
|
||||||
@@ -60,13 +60,13 @@ static void PAConvectionSetup2D(const int NQ,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA Convection Assemble 3D kernel
|
// PA Convection Assemble 3D kernel
|
||||||
static void PAConvectionSetup3D(const int NQ,
|
void PAConvectionSetup3D(const int NQ,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
const Vector &vel,
|
const Vector &vel,
|
||||||
const double alpha,
|
const double alpha,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
constexpr int DIM = 3;
|
constexpr int DIM = 3;
|
||||||
constexpr int SDIM = DIM;
|
constexpr int SDIM = DIM;
|
||||||
@@ -135,7 +135,7 @@ static void PAConvectionSetup(const int dim,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA Convection Apply 2D kernel
|
// PA Convection Apply 2D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0> static
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
void PAConvectionApply2D(const int ne,
|
void PAConvectionApply2D(const int ne,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
@@ -254,7 +254,7 @@ void PAConvectionApply2D(const int ne,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// Optimized PA Convection Apply 2D kernel
|
// Optimized PA Convection Apply 2D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0> static
|
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0>
|
||||||
void SmemPAConvectionApply2D(const int ne,
|
void SmemPAConvectionApply2D(const int ne,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
@@ -382,7 +382,7 @@ void SmemPAConvectionApply2D(const int ne,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA Convection Apply 3D kernel
|
// PA Convection Apply 3D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0> static
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
void PAConvectionApply3D(const int ne,
|
void PAConvectionApply3D(const int ne,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
@@ -563,7 +563,7 @@ void PAConvectionApply3D(const int ne,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// Optimized PA Convection Apply 3D kernel
|
// Optimized PA Convection Apply 3D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0> static
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
void SmemPAConvectionApply3D(const int ne,
|
void SmemPAConvectionApply3D(const int ne,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
|
|||||||
@@ -16,12 +16,12 @@
|
|||||||
namespace mfem
|
namespace mfem
|
||||||
{
|
{
|
||||||
|
|
||||||
static void EADGTraceAssemble1DInt(const int NF,
|
void EADGTraceAssemble1DInt(const int NF,
|
||||||
const Array<double> &basis,
|
const Array<double> &basis,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata_int,
|
Vector &eadata_int,
|
||||||
Vector &eadata_ext,
|
Vector &eadata_ext,
|
||||||
const bool add)
|
const bool add)
|
||||||
{
|
{
|
||||||
auto D = Reshape(padata.Read(), 2, 2, NF);
|
auto D = Reshape(padata.Read(), 2, 2, NF);
|
||||||
auto A_int = Reshape(eadata_int.ReadWrite(), 2, NF);
|
auto A_int = Reshape(eadata_int.ReadWrite(), 2, NF);
|
||||||
@@ -50,11 +50,11 @@ static void EADGTraceAssemble1DInt(const int NF,
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void EADGTraceAssemble1DBdr(const int NF,
|
void EADGTraceAssemble1DBdr(const int NF,
|
||||||
const Array<double> &basis,
|
const Array<double> &basis,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata_bdr,
|
Vector &eadata_bdr,
|
||||||
const bool add)
|
const bool add)
|
||||||
{
|
{
|
||||||
auto D = Reshape(padata.Read(), 2, 2, NF);
|
auto D = Reshape(padata.Read(), 2, 2, NF);
|
||||||
auto A_bdr = Reshape(eadata_bdr.ReadWrite(), NF);
|
auto A_bdr = Reshape(eadata_bdr.ReadWrite(), NF);
|
||||||
@@ -72,14 +72,14 @@ static void EADGTraceAssemble1DBdr(const int NF,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EADGTraceAssemble2DInt(const int NF,
|
void EADGTraceAssemble2DInt(const int NF,
|
||||||
const Array<double> &basis,
|
const Array<double> &basis,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata_int,
|
Vector &eadata_int,
|
||||||
Vector &eadata_ext,
|
Vector &eadata_ext,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -128,13 +128,13 @@ static void EADGTraceAssemble2DInt(const int NF,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EADGTraceAssemble2DBdr(const int NF,
|
void EADGTraceAssemble2DBdr(const int NF,
|
||||||
const Array<double> &basis,
|
const Array<double> &basis,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata_bdr,
|
Vector &eadata_bdr,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -170,14 +170,14 @@ static void EADGTraceAssemble2DBdr(const int NF,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EADGTraceAssemble3DInt(const int NF,
|
void EADGTraceAssemble3DInt(const int NF,
|
||||||
const Array<double> &basis,
|
const Array<double> &basis,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata_int,
|
Vector &eadata_int,
|
||||||
Vector &eadata_ext,
|
Vector &eadata_ext,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -268,13 +268,13 @@ static void EADGTraceAssemble3DInt(const int NF,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EADGTraceAssemble3DBdr(const int NF,
|
void EADGTraceAssemble3DBdr(const int NF,
|
||||||
const Array<double> &basis,
|
const Array<double> &basis,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata_bdr,
|
Vector &eadata_bdr,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
|
|||||||
@@ -19,16 +19,16 @@ using namespace std;
|
|||||||
namespace mfem
|
namespace mfem
|
||||||
{
|
{
|
||||||
// PA DG Trace Integrator
|
// PA DG Trace Integrator
|
||||||
static void PADGTraceSetup2D(const int Q1D,
|
void PADGTraceSetup2D(const int Q1D,
|
||||||
const int NF,
|
const int NF,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &det,
|
const Vector &det,
|
||||||
const Vector &nor,
|
const Vector &nor,
|
||||||
const Vector &rho,
|
const Vector &rho,
|
||||||
const Vector &vel,
|
const Vector &vel,
|
||||||
const double alpha,
|
const double alpha,
|
||||||
const double beta,
|
const double beta,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int VDIM = 2;
|
const int VDIM = 2;
|
||||||
|
|
||||||
@@ -61,16 +61,16 @@ static void PADGTraceSetup2D(const int Q1D,
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PADGTraceSetup3D(const int Q1D,
|
void PADGTraceSetup3D(const int Q1D,
|
||||||
const int NF,
|
const int NF,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &det,
|
const Vector &det,
|
||||||
const Vector &nor,
|
const Vector &nor,
|
||||||
const Vector &rho,
|
const Vector &rho,
|
||||||
const Vector &vel,
|
const Vector &vel,
|
||||||
const double alpha,
|
const double alpha,
|
||||||
const double beta,
|
const double beta,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int VDIM = 3;
|
const int VDIM = 3;
|
||||||
|
|
||||||
@@ -301,7 +301,7 @@ void DGTraceIntegrator::AssemblePABoundaryFaces(const FiniteElementSpace& fes)
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA DGTrace Apply 2D kernel for Gauss-Lobatto/Bernstein
|
// PA DGTrace Apply 2D kernel for Gauss-Lobatto/Bernstein
|
||||||
template<int T_D1D = 0, int T_Q1D = 0> static
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
void PADGTraceApply2D(const int NF,
|
void PADGTraceApply2D(const int NF,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
@@ -392,7 +392,7 @@ void PADGTraceApply2D(const int NF,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA DGTrace Apply 3D kernel for Gauss-Lobatto/Bernstein
|
// PA DGTrace Apply 3D kernel for Gauss-Lobatto/Bernstein
|
||||||
template<int T_D1D = 0, int T_Q1D = 0> static
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
void PADGTraceApply3D(const int NF,
|
void PADGTraceApply3D(const int NF,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
@@ -537,7 +537,7 @@ void PADGTraceApply3D(const int NF,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// Optimized PA DGTrace Apply 3D kernel for Gauss-Lobatto/Bernstein
|
// Optimized PA DGTrace Apply 3D kernel for Gauss-Lobatto/Bernstein
|
||||||
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0> static
|
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0>
|
||||||
void SmemPADGTraceApply3D(const int NF,
|
void SmemPADGTraceApply3D(const int NF,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
@@ -701,7 +701,7 @@ static void PADGTraceApply(const int dim,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA DGTrace Apply 2D kernel for Gauss-Lobatto/Bernstein
|
// PA DGTrace Apply 2D kernel for Gauss-Lobatto/Bernstein
|
||||||
template<int T_D1D = 0, int T_Q1D = 0> static
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
void PADGTraceApplyTranspose2D(const int NF,
|
void PADGTraceApplyTranspose2D(const int NF,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
@@ -797,7 +797,7 @@ void PADGTraceApplyTranspose2D(const int NF,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA DGTrace Apply Transpose 3D kernel for Gauss-Lobatto/Bernstein
|
// PA DGTrace Apply Transpose 3D kernel for Gauss-Lobatto/Bernstein
|
||||||
template<int T_D1D = 0, int T_Q1D = 0> static
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
void PADGTraceApplyTranspose3D(const int NF,
|
void PADGTraceApplyTranspose3D(const int NF,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
@@ -953,7 +953,7 @@ void PADGTraceApplyTranspose3D(const int NF,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// Optimized PA DGTrace Apply Transpose 3D kernel for Gauss-Lobatto/Bernstein
|
// Optimized PA DGTrace Apply Transpose 3D kernel for Gauss-Lobatto/Bernstein
|
||||||
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0> static
|
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0>
|
||||||
void SmemPADGTraceApplyTranspose3D(const int NF,
|
void SmemPADGTraceApplyTranspose3D(const int NF,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
|
|||||||
@@ -17,14 +17,14 @@ namespace mfem
|
|||||||
{
|
{
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EADiffusionAssemble1D(const int NE,
|
void EADiffusionAssemble1D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata,
|
Vector &eadata,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -68,14 +68,14 @@ static void EADiffusionAssemble1D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EADiffusionAssemble2D(const int NE,
|
void EADiffusionAssemble2D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata,
|
Vector &eadata,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -145,14 +145,14 @@ static void EADiffusionAssemble2D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EADiffusionAssemble3D(const int NE,
|
void EADiffusionAssemble3D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata,
|
Vector &eadata,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
|
|||||||
@@ -496,14 +496,14 @@ void DiffusionIntegrator::AssemblePA(const FiniteElementSpace &fes)
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PADiffusionDiagonal2D(const int NE,
|
void PADiffusionDiagonal2D(const int NE,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Vector &d,
|
const Vector &d,
|
||||||
Vector &y,
|
Vector &y,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -562,14 +562,14 @@ static void PADiffusionDiagonal2D(const int NE,
|
|||||||
|
|
||||||
// Shared memory PA Diffusion Diagonal 2D kernel
|
// Shared memory PA Diffusion Diagonal 2D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0>
|
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0>
|
||||||
static void SmemPADiffusionDiagonal2D(const int NE,
|
void SmemPADiffusionDiagonal2D(const int NE,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &g_,
|
const Array<double> &g_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -656,14 +656,14 @@ static void SmemPADiffusionDiagonal2D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PADiffusionDiagonal3D(const int NE,
|
void PADiffusionDiagonal3D(const int NE,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Vector &d,
|
const Vector &d,
|
||||||
Vector &y,
|
Vector &y,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
constexpr int DIM = 3;
|
constexpr int DIM = 3;
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
@@ -757,14 +757,14 @@ static void PADiffusionDiagonal3D(const int NE,
|
|||||||
|
|
||||||
// Shared memory PA Diffusion Diagonal 3D kernel
|
// Shared memory PA Diffusion Diagonal 3D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void SmemPADiffusionDiagonal3D(const int NE,
|
void SmemPADiffusionDiagonal3D(const int NE,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &g_,
|
const Array<double> &g_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
constexpr int DIM = 3;
|
constexpr int DIM = 3;
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
@@ -1034,17 +1034,17 @@ static void OccaPADiffusionApply3D(const int D1D,
|
|||||||
|
|
||||||
// PA Diffusion Apply 2D kernel
|
// PA Diffusion Apply 2D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PADiffusionApply2D(const int NE,
|
void PADiffusionApply2D(const int NE,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &g_,
|
const Array<double> &g_,
|
||||||
const Array<double> &bt_,
|
const Array<double> &bt_,
|
||||||
const Array<double> >_,
|
const Array<double> >_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -1156,15 +1156,15 @@ static void PADiffusionApply2D(const int NE,
|
|||||||
|
|
||||||
// Shared memory PA Diffusion Apply 2D kernel
|
// Shared memory PA Diffusion Apply 2D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0>
|
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0>
|
||||||
static void SmemPADiffusionApply2D(const int NE,
|
void SmemPADiffusionApply2D(const int NE,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &g_,
|
const Array<double> &g_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -1314,16 +1314,16 @@ static void SmemPADiffusionApply2D(const int NE,
|
|||||||
|
|
||||||
// PA Diffusion Apply 3D kernel
|
// PA Diffusion Apply 3D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PADiffusionApply3D(const int NE,
|
void PADiffusionApply3D(const int NE,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
const Array<double> >,
|
const Array<double> >,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
int d1d = 0, int q1d = 0)
|
int d1d = 0, int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -1533,15 +1533,15 @@ static MFEM_HOST_DEVICE inline double sign(const int q, const int d)
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void SmemPADiffusionApply3D(const int NE,
|
void SmemPADiffusionApply3D(const int NE,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &g_,
|
const Array<double> &g_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
|
|||||||
@@ -21,12 +21,12 @@ namespace mfem
|
|||||||
// PA Divergence Integrator
|
// PA Divergence Integrator
|
||||||
|
|
||||||
// PA Divergence Assemble 2D kernel
|
// PA Divergence Assemble 2D kernel
|
||||||
static void PADivergenceSetup2D(const int Q1D,
|
void PADivergenceSetup2D(const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
const double COEFF,
|
const double COEFF,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D;
|
const int NQ = Q1D*Q1D;
|
||||||
auto W = w.Read();
|
auto W = w.Read();
|
||||||
@@ -51,12 +51,12 @@ static void PADivergenceSetup2D(const int Q1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA Divergence Assemble 3D kernel
|
// PA Divergence Assemble 3D kernel
|
||||||
static void PADivergenceSetup3D(const int Q1D,
|
void PADivergenceSetup3D(const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
const double COEFF,
|
const double COEFF,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D*Q1D;
|
const int NQ = Q1D*Q1D*Q1D;
|
||||||
auto W = w.Read();
|
auto W = w.Read();
|
||||||
@@ -160,16 +160,16 @@ void VectorDivergenceIntegrator::AssemblePA(const FiniteElementSpace &trial_fes,
|
|||||||
|
|
||||||
// PA Divergence Apply 2D kernel
|
// PA Divergence Apply 2D kernel
|
||||||
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
||||||
static void PADivergenceApply2D(const int NE,
|
void PADivergenceApply2D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int tr_d1d = 0,
|
const int tr_d1d = 0,
|
||||||
const int te_d1d = 0,
|
const int te_d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
||||||
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
||||||
@@ -281,16 +281,16 @@ static void PADivergenceApply2D(const int NE,
|
|||||||
// Shared memory PA Divergence Apply 2D kernel
|
// Shared memory PA Divergence Apply 2D kernel
|
||||||
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0,
|
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0,
|
||||||
const int T_NBZ = 0>
|
const int T_NBZ = 0>
|
||||||
static void SmemPADivergenceApply2D(const int NE,
|
void SmemPADivergenceApply2D(const int NE,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &g_,
|
const Array<double> &g_,
|
||||||
const Array<double> &bt_,
|
const Array<double> &bt_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int tr_d1d = 0,
|
const int tr_d1d = 0,
|
||||||
const int te_d1d = 0,
|
const int te_d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
// TODO
|
// TODO
|
||||||
MFEM_ASSERT(false, "SHARED MEM NOT PROGRAMMED YET");
|
MFEM_ASSERT(false, "SHARED MEM NOT PROGRAMMED YET");
|
||||||
@@ -298,16 +298,16 @@ static void SmemPADivergenceApply2D(const int NE,
|
|||||||
|
|
||||||
// PA Divergence Apply 2D kernel transpose
|
// PA Divergence Apply 2D kernel transpose
|
||||||
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
||||||
static void PADivergenceApplyTranspose2D(const int NE,
|
void PADivergenceApplyTranspose2D(const int NE,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
const Array<double> >,
|
const Array<double> >,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int tr_d1d = 0,
|
const int tr_d1d = 0,
|
||||||
const int te_d1d = 0,
|
const int te_d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
||||||
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
||||||
@@ -414,16 +414,16 @@ static void PADivergenceApplyTranspose2D(const int NE,
|
|||||||
|
|
||||||
// PA Vector Divergence Apply 3D kernel
|
// PA Vector Divergence Apply 3D kernel
|
||||||
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
||||||
static void PADivergenceApply3D(const int NE,
|
void PADivergenceApply3D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
int tr_d1d = 0,
|
int tr_d1d = 0,
|
||||||
int te_d1d = 0,
|
int te_d1d = 0,
|
||||||
int q1d = 0)
|
int q1d = 0)
|
||||||
{
|
{
|
||||||
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
||||||
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
||||||
@@ -597,16 +597,16 @@ static void PADivergenceApply3D(const int NE,
|
|||||||
|
|
||||||
// PA Vector Divergence Apply 3D kernel
|
// PA Vector Divergence Apply 3D kernel
|
||||||
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
||||||
static void PADivergenceApplyTranspose3D(const int NE,
|
void PADivergenceApplyTranspose3D(const int NE,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
const Array<double> >,
|
const Array<double> >,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
int tr_d1d = 0,
|
int tr_d1d = 0,
|
||||||
int te_d1d = 0,
|
int te_d1d = 0,
|
||||||
int q1d = 0)
|
int q1d = 0)
|
||||||
{
|
{
|
||||||
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
||||||
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
||||||
@@ -775,16 +775,16 @@ static void PADivergenceApplyTranspose3D(const int NE,
|
|||||||
|
|
||||||
// Shared memory PA Vector Divergence Apply 3D kernel
|
// Shared memory PA Vector Divergence Apply 3D kernel
|
||||||
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
||||||
static void SmemPADivergenceApply3D(const int NE,
|
void SmemPADivergenceApply3D(const int NE,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &g_,
|
const Array<double> &g_,
|
||||||
const Array<double> &bt_,
|
const Array<double> &bt_,
|
||||||
const Vector &q_,
|
const Vector &q_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int tr_d1d = 0,
|
const int tr_d1d = 0,
|
||||||
const int te_d1d = 0,
|
const int te_d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
||||||
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
||||||
|
|||||||
+42
-42
@@ -70,12 +70,12 @@ namespace mfem
|
|||||||
the \b MFEM_SHARED keyword for local arrays. */
|
the \b MFEM_SHARED keyword for local arrays. */
|
||||||
|
|
||||||
// PA Gradient Assemble 2D kernel
|
// PA Gradient Assemble 2D kernel
|
||||||
static void PAGradientSetup2D(const int Q1D,
|
void PAGradientSetup2D(const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
const Vector &c,
|
const Vector &c,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D;
|
const int NQ = Q1D*Q1D;
|
||||||
auto W = w.Read();
|
auto W = w.Read();
|
||||||
@@ -105,12 +105,12 @@ static void PAGradientSetup2D(const int Q1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA Gradient Assemble 3D kernel
|
// PA Gradient Assemble 3D kernel
|
||||||
static void PAGradientSetup3D(const int Q1D,
|
void PAGradientSetup3D(const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
const Vector &c,
|
const Vector &c,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D*Q1D;
|
const int NQ = Q1D*Q1D*Q1D;
|
||||||
auto W = w.Read();
|
auto W = w.Read();
|
||||||
@@ -254,16 +254,16 @@ void GradientIntegrator::AssemblePA(const FiniteElementSpace &trial_fes,
|
|||||||
|
|
||||||
// PA Gradient Apply 2D kernel
|
// PA Gradient Apply 2D kernel
|
||||||
template<int T_TR_D1D = 0, int T_TE_D1D = 0, int T_Q1D = 0>
|
template<int T_TR_D1D = 0, int T_TE_D1D = 0, int T_Q1D = 0>
|
||||||
static void PAGradientApply2D(const int NE,
|
void PAGradientApply2D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int tr_d1d = 0,
|
const int tr_d1d = 0,
|
||||||
const int te_d1d = 0,
|
const int te_d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
||||||
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
||||||
@@ -384,16 +384,16 @@ static void PAGradientApplyTranspose2D(const int NE,
|
|||||||
|
|
||||||
// PA Gradient Apply 3D kernel
|
// PA Gradient Apply 3D kernel
|
||||||
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
||||||
static void PAGradientApply3D(const int NE,
|
void PAGradientApply3D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
int tr_d1d = 0,
|
int tr_d1d = 0,
|
||||||
int te_d1d = 0,
|
int te_d1d = 0,
|
||||||
int q1d = 0)
|
int q1d = 0)
|
||||||
{
|
{
|
||||||
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
||||||
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
||||||
@@ -579,16 +579,16 @@ static void PAGradientApplyTranspose3D(const int NE,
|
|||||||
|
|
||||||
// Shared memory PA Gradient Apply 3D kernel
|
// Shared memory PA Gradient Apply 3D kernel
|
||||||
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
template<const int T_TR_D1D = 0, const int T_TE_D1D = 0, const int T_Q1D = 0>
|
||||||
static void SmemPAGradientApply3D(const int NE,
|
void SmemPAGradientApply3D(const int NE,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &g_,
|
const Array<double> &g_,
|
||||||
const Array<double> &bt_,
|
const Array<double> &bt_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int tr_d1d = 0,
|
const int tr_d1d = 0,
|
||||||
const int te_d1d = 0,
|
const int te_d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
const int TR_D1D = T_TR_D1D ? T_TR_D1D : tr_d1d;
|
||||||
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
const int TE_D1D = T_TE_D1D ? T_TE_D1D : te_d1d;
|
||||||
|
|||||||
+194
-194
@@ -791,12 +791,12 @@ void SmemPAHcurlMassApply3D(const int D1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA H(curl) curl-curl assemble 2D kernel
|
// PA H(curl) curl-curl assemble 2D kernel
|
||||||
static void PACurlCurlSetup2D(const int Q1D,
|
void PACurlCurlSetup2D(const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
Vector &coeff,
|
Vector &coeff,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D;
|
const int NQ = Q1D*Q1D;
|
||||||
auto W = w.Read();
|
auto W = w.Read();
|
||||||
@@ -818,13 +818,13 @@ static void PACurlCurlSetup2D(const int Q1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA H(curl) curl-curl assemble 3D kernel
|
// PA H(curl) curl-curl assemble 3D kernel
|
||||||
static void PACurlCurlSetup3D(const int Q1D,
|
void PACurlCurlSetup3D(const int Q1D,
|
||||||
const int coeffDim,
|
const int coeffDim,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
Vector &coeff,
|
Vector &coeff,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D*Q1D;
|
const int NQ = Q1D*Q1D*Q1D;
|
||||||
const bool symmetric = (coeffDim != 9);
|
const bool symmetric = (coeffDim != 9);
|
||||||
@@ -1045,16 +1045,16 @@ void CurlCurlIntegrator::AssemblePA(const FiniteElementSpace &fes)
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PACurlCurlApply2D(const int D1D,
|
void PACurlCurlApply2D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &bo,
|
const Array<double> &bo,
|
||||||
const Array<double> &bot,
|
const Array<double> &bot,
|
||||||
const Array<double> &gc,
|
const Array<double> &gc,
|
||||||
const Array<double> &gct,
|
const Array<double> &gct,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x,
|
const Vector &x,
|
||||||
Vector &y)
|
Vector &y)
|
||||||
{
|
{
|
||||||
constexpr static int VDIM = 2;
|
constexpr static int VDIM = 2;
|
||||||
constexpr static int MAX_D1D = HCURL_MAX_D1D;
|
constexpr static int MAX_D1D = HCURL_MAX_D1D;
|
||||||
@@ -1166,19 +1166,19 @@ static void PACurlCurlApply2D(const int D1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
||||||
static void PACurlCurlApply3D(const int D1D,
|
void PACurlCurlApply3D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &bo,
|
const Array<double> &bo,
|
||||||
const Array<double> &bc,
|
const Array<double> &bc,
|
||||||
const Array<double> &bot,
|
const Array<double> &bot,
|
||||||
const Array<double> &bct,
|
const Array<double> &bct,
|
||||||
const Array<double> &gc,
|
const Array<double> &gc,
|
||||||
const Array<double> &gct,
|
const Array<double> &gct,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x,
|
const Vector &x,
|
||||||
Vector &y)
|
Vector &y)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
||||||
@@ -1677,19 +1677,19 @@ static void PACurlCurlApply3D(const int D1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
||||||
static void SmemPACurlCurlApply3D(const int D1D,
|
void SmemPACurlCurlApply3D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &bo,
|
const Array<double> &bo,
|
||||||
const Array<double> &bc,
|
const Array<double> &bc,
|
||||||
const Array<double> &bot,
|
const Array<double> &bot,
|
||||||
const Array<double> &bct,
|
const Array<double> &bct,
|
||||||
const Array<double> &gc,
|
const Array<double> &gc,
|
||||||
const Array<double> &gct,
|
const Array<double> &gct,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x,
|
const Vector &x,
|
||||||
Vector &y)
|
Vector &y)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
||||||
@@ -2032,13 +2032,13 @@ void CurlCurlIntegrator::AddMultPA(const Vector &x, Vector &y) const
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PACurlCurlAssembleDiagonal2D(const int D1D,
|
void PACurlCurlAssembleDiagonal2D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &bo,
|
const Array<double> &bo,
|
||||||
const Array<double> &gc,
|
const Array<double> &gc,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
Vector &diag)
|
Vector &diag)
|
||||||
{
|
{
|
||||||
constexpr static int VDIM = 2;
|
constexpr static int VDIM = 2;
|
||||||
constexpr static int MAX_Q1D = HCURL_MAX_Q1D;
|
constexpr static int MAX_Q1D = HCURL_MAX_Q1D;
|
||||||
@@ -2087,16 +2087,16 @@ static void PACurlCurlAssembleDiagonal2D(const int D1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
||||||
static void PACurlCurlAssembleDiagonal3D(const int D1D,
|
void PACurlCurlAssembleDiagonal3D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &bo,
|
const Array<double> &bo,
|
||||||
const Array<double> &bc,
|
const Array<double> &bc,
|
||||||
const Array<double> &go,
|
const Array<double> &go,
|
||||||
const Array<double> &gc,
|
const Array<double> &gc,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
Vector &diag)
|
Vector &diag)
|
||||||
{
|
{
|
||||||
constexpr static int VDIM = 3;
|
constexpr static int VDIM = 3;
|
||||||
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
||||||
@@ -2273,16 +2273,16 @@ static void PACurlCurlAssembleDiagonal3D(const int D1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
||||||
static void SmemPACurlCurlAssembleDiagonal3D(const int D1D,
|
void SmemPACurlCurlAssembleDiagonal3D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const bool symmetric,
|
const bool symmetric,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &bo,
|
const Array<double> &bo,
|
||||||
const Array<double> &bc,
|
const Array<double> &bc,
|
||||||
const Array<double> &go,
|
const Array<double> &go,
|
||||||
const Array<double> &gc,
|
const Array<double> &gc,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
Vector &diag)
|
Vector &diag)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
||||||
@@ -2955,18 +2955,18 @@ void MixedVectorCurlIntegrator::AssemblePA(const FiniteElementSpace &trial_fes,
|
|||||||
// Apply to x corresponding to DOF's in H(curl) (trial), whose curl is
|
// Apply to x corresponding to DOF's in H(curl) (trial), whose curl is
|
||||||
// integrated against H(curl) test functions corresponding to y.
|
// integrated against H(curl) test functions corresponding to y.
|
||||||
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
||||||
static void PAHcurlL2Apply3D(const int D1D,
|
void PAHcurlL2Apply3D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int coeffDim,
|
const int coeffDim,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &bo,
|
const Array<double> &bo,
|
||||||
const Array<double> &bc,
|
const Array<double> &bc,
|
||||||
const Array<double> &bot,
|
const Array<double> &bot,
|
||||||
const Array<double> &bct,
|
const Array<double> &bct,
|
||||||
const Array<double> &gc,
|
const Array<double> &gc,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x,
|
const Vector &x,
|
||||||
Vector &y)
|
Vector &y)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
||||||
@@ -3297,16 +3297,16 @@ static void PAHcurlL2Apply3D(const int D1D,
|
|||||||
// Apply to x corresponding to DOF's in H(curl) (trial), whose curl is
|
// Apply to x corresponding to DOF's in H(curl) (trial), whose curl is
|
||||||
// integrated against H(curl) test functions corresponding to y.
|
// integrated against H(curl) test functions corresponding to y.
|
||||||
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
||||||
static void SmemPAHcurlL2Apply3D(const int D1D,
|
void SmemPAHcurlL2Apply3D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int coeffDim,
|
const int coeffDim,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &bo,
|
const Array<double> &bo,
|
||||||
const Array<double> &bc,
|
const Array<double> &bc,
|
||||||
const Array<double> &gc,
|
const Array<double> &gc,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x,
|
const Vector &x,
|
||||||
Vector &y)
|
Vector &y)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
||||||
@@ -3585,18 +3585,18 @@ static void SmemPAHcurlL2Apply3D(const int D1D,
|
|||||||
// Apply to x corresponding to DOF's in H(curl) (trial), whose curl is
|
// Apply to x corresponding to DOF's in H(curl) (trial), whose curl is
|
||||||
// integrated against H(div) test functions corresponding to y.
|
// integrated against H(div) test functions corresponding to y.
|
||||||
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
||||||
static void PAHcurlHdivApply3D(const int D1D,
|
void PAHcurlHdivApply3D(const int D1D,
|
||||||
const int D1Dtest,
|
const int D1Dtest,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &bo,
|
const Array<double> &bo,
|
||||||
const Array<double> &bc,
|
const Array<double> &bc,
|
||||||
const Array<double> &bot,
|
const Array<double> &bot,
|
||||||
const Array<double> &bct,
|
const Array<double> &bct,
|
||||||
const Array<double> &gc,
|
const Array<double> &gc,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x,
|
const Vector &x,
|
||||||
Vector &y)
|
Vector &y)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
||||||
@@ -4071,18 +4071,18 @@ void MixedVectorWeakCurlIntegrator::AssemblePA(const FiniteElementSpace
|
|||||||
// Apply to x corresponding to DOF's in H(curl) (trial), integrated against curl
|
// Apply to x corresponding to DOF's in H(curl) (trial), integrated against curl
|
||||||
// of H(curl) test functions corresponding to y.
|
// of H(curl) test functions corresponding to y.
|
||||||
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
||||||
static void PAHcurlL2Apply3DTranspose(const int D1D,
|
void PAHcurlL2Apply3DTranspose(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int coeffDim,
|
const int coeffDim,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &bo,
|
const Array<double> &bo,
|
||||||
const Array<double> &bc,
|
const Array<double> &bc,
|
||||||
const Array<double> &bot,
|
const Array<double> &bot,
|
||||||
const Array<double> &bct,
|
const Array<double> &bct,
|
||||||
const Array<double> &gct,
|
const Array<double> &gct,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x,
|
const Vector &x,
|
||||||
Vector &y)
|
Vector &y)
|
||||||
{
|
{
|
||||||
// See PAHcurlL2Apply3D for comments.
|
// See PAHcurlL2Apply3D for comments.
|
||||||
|
|
||||||
@@ -4413,16 +4413,16 @@ static void PAHcurlL2Apply3DTranspose(const int D1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
template<int MAX_D1D = HCURL_MAX_D1D, int MAX_Q1D = HCURL_MAX_Q1D>
|
||||||
static void SmemPAHcurlL2Apply3DTranspose(const int D1D,
|
void SmemPAHcurlL2Apply3DTranspose(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int coeffDim,
|
const int coeffDim,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &bo,
|
const Array<double> &bo,
|
||||||
const Array<double> &bc,
|
const Array<double> &bc,
|
||||||
const Array<double> &gc,
|
const Array<double> &gc,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x,
|
const Vector &x,
|
||||||
Vector &y)
|
Vector &y)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
MFEM_VERIFY(D1D <= MAX_D1D, "Error: D1D > MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
MFEM_VERIFY(Q1D <= MAX_Q1D, "Error: Q1D > MAX_Q1D");
|
||||||
@@ -4675,13 +4675,13 @@ void MixedVectorWeakCurlIntegrator::AddMultPA(const Vector &x, Vector &y) const
|
|||||||
// Apply to x corresponding to DOFs in H^1 (domain) the (topological) gradient
|
// Apply to x corresponding to DOFs in H^1 (domain) the (topological) gradient
|
||||||
// to get a dof in H(curl) (range). You can think of the range as the "test" space
|
// to get a dof in H(curl) (range). You can think of the range as the "test" space
|
||||||
// and the domain as the "trial" space, but there's no integration.
|
// and the domain as the "trial" space, but there's no integration.
|
||||||
static void PAHcurlApplyGradient2D(const int c_dofs1D,
|
void PAHcurlApplyGradient2D(const int c_dofs1D,
|
||||||
const int o_dofs1D,
|
const int o_dofs1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &B_,
|
const Array<double> &B_,
|
||||||
const Array<double> &G_,
|
const Array<double> &G_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
auto B = Reshape(B_.Read(), c_dofs1D, c_dofs1D);
|
auto B = Reshape(B_.Read(), c_dofs1D, c_dofs1D);
|
||||||
auto G = Reshape(G_.Read(), o_dofs1D, c_dofs1D);
|
auto G = Reshape(G_.Read(), o_dofs1D, c_dofs1D);
|
||||||
@@ -4753,12 +4753,12 @@ static void PAHcurlApplyGradient2D(const int c_dofs1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// Specialization of PAHcurlApplyGradient2D to the case where B is identity
|
// Specialization of PAHcurlApplyGradient2D to the case where B is identity
|
||||||
static void PAHcurlApplyGradient2DBId(const int c_dofs1D,
|
void PAHcurlApplyGradient2DBId(const int c_dofs1D,
|
||||||
const int o_dofs1D,
|
const int o_dofs1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &G_,
|
const Array<double> &G_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
auto G = Reshape(G_.Read(), o_dofs1D, c_dofs1D);
|
auto G = Reshape(G_.Read(), o_dofs1D, c_dofs1D);
|
||||||
|
|
||||||
@@ -4822,7 +4822,7 @@ static void PAHcurlApplyGradient2DBId(const int c_dofs1D,
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PAHcurlApplyGradientTranspose2D(
|
void PAHcurlApplyGradientTranspose2D(
|
||||||
const int c_dofs1D, const int o_dofs1D, const int NE,
|
const int c_dofs1D, const int o_dofs1D, const int NE,
|
||||||
const Array<double> &B_, const Array<double> &G_,
|
const Array<double> &B_, const Array<double> &G_,
|
||||||
const Vector &x_, Vector &y_)
|
const Vector &x_, Vector &y_)
|
||||||
@@ -4898,7 +4898,7 @@ static void PAHcurlApplyGradientTranspose2D(
|
|||||||
|
|
||||||
// Specialization of PAHcurlApplyGradientTranspose2D to the case where
|
// Specialization of PAHcurlApplyGradientTranspose2D to the case where
|
||||||
// B is identity
|
// B is identity
|
||||||
static void PAHcurlApplyGradientTranspose2DBId(
|
void PAHcurlApplyGradientTranspose2DBId(
|
||||||
const int c_dofs1D, const int o_dofs1D, const int NE,
|
const int c_dofs1D, const int o_dofs1D, const int NE,
|
||||||
const Array<double> &G_,
|
const Array<double> &G_,
|
||||||
const Vector &x_, Vector &y_)
|
const Vector &x_, Vector &y_)
|
||||||
@@ -4965,13 +4965,13 @@ static void PAHcurlApplyGradientTranspose2DBId(
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PAHcurlApplyGradient3D(const int c_dofs1D,
|
void PAHcurlApplyGradient3D(const int c_dofs1D,
|
||||||
const int o_dofs1D,
|
const int o_dofs1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &B_,
|
const Array<double> &B_,
|
||||||
const Array<double> &G_,
|
const Array<double> &G_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
auto B = Reshape(B_.Read(), c_dofs1D, c_dofs1D);
|
auto B = Reshape(B_.Read(), c_dofs1D, c_dofs1D);
|
||||||
auto G = Reshape(G_.Read(), o_dofs1D, c_dofs1D);
|
auto G = Reshape(G_.Read(), o_dofs1D, c_dofs1D);
|
||||||
@@ -5154,12 +5154,12 @@ static void PAHcurlApplyGradient3D(const int c_dofs1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// Specialization of PAHcurlApplyGradient3D to the case where
|
// Specialization of PAHcurlApplyGradient3D to the case where
|
||||||
static void PAHcurlApplyGradient3DBId(const int c_dofs1D,
|
void PAHcurlApplyGradient3DBId(const int c_dofs1D,
|
||||||
const int o_dofs1D,
|
const int o_dofs1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &G_,
|
const Array<double> &G_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
auto G = Reshape(G_.Read(), o_dofs1D, c_dofs1D);
|
auto G = Reshape(G_.Read(), o_dofs1D, c_dofs1D);
|
||||||
|
|
||||||
@@ -5322,7 +5322,7 @@ static void PAHcurlApplyGradient3DBId(const int c_dofs1D,
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PAHcurlApplyGradientTranspose3D(
|
void PAHcurlApplyGradientTranspose3D(
|
||||||
const int c_dofs1D, const int o_dofs1D, const int NE,
|
const int c_dofs1D, const int o_dofs1D, const int NE,
|
||||||
const Array<double> &B_, const Array<double> &G_,
|
const Array<double> &B_, const Array<double> &G_,
|
||||||
const Vector &x_, Vector &y_)
|
const Vector &x_, Vector &y_)
|
||||||
@@ -5507,7 +5507,7 @@ static void PAHcurlApplyGradientTranspose3D(
|
|||||||
}
|
}
|
||||||
|
|
||||||
// Specialization of PAHcurlApplyGradientTranspose3D to the case where
|
// Specialization of PAHcurlApplyGradientTranspose3D to the case where
|
||||||
static void PAHcurlApplyGradientTranspose3DBId(
|
void PAHcurlApplyGradientTranspose3DBId(
|
||||||
const int c_dofs1D, const int o_dofs1D, const int NE,
|
const int c_dofs1D, const int o_dofs1D, const int NE,
|
||||||
const Array<double> &G_,
|
const Array<double> &G_,
|
||||||
const Vector &x_, Vector &y_)
|
const Vector &x_, Vector &y_)
|
||||||
@@ -5789,14 +5789,14 @@ void GradientInterpolator::AddMultTransposePA(const Vector &x, Vector &y) const
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PAHcurlVecH1IdentityApply3D(const int c_dofs1D,
|
void PAHcurlVecH1IdentityApply3D(const int c_dofs1D,
|
||||||
const int o_dofs1D,
|
const int o_dofs1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &Bclosed,
|
const Array<double> &Bclosed,
|
||||||
const Array<double> &Bopen,
|
const Array<double> &Bopen,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
auto Bc = Reshape(Bclosed.Read(), c_dofs1D, c_dofs1D);
|
auto Bc = Reshape(Bclosed.Read(), c_dofs1D, c_dofs1D);
|
||||||
auto Bo = Reshape(Bopen.Read(), o_dofs1D, c_dofs1D);
|
auto Bo = Reshape(Bopen.Read(), o_dofs1D, c_dofs1D);
|
||||||
@@ -6002,14 +6002,14 @@ static void PAHcurlVecH1IdentityApply3D(const int c_dofs1D,
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PAHcurlVecH1IdentityApplyTranspose3D(const int c_dofs1D,
|
void PAHcurlVecH1IdentityApplyTranspose3D(const int c_dofs1D,
|
||||||
const int o_dofs1D,
|
const int o_dofs1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &Bclosed,
|
const Array<double> &Bclosed,
|
||||||
const Array<double> &Bopen,
|
const Array<double> &Bopen,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
auto Bc = Reshape(Bclosed.Read(), c_dofs1D, c_dofs1D);
|
auto Bc = Reshape(Bclosed.Read(), c_dofs1D, c_dofs1D);
|
||||||
auto Bo = Reshape(Bopen.Read(), o_dofs1D, c_dofs1D);
|
auto Bo = Reshape(Bopen.Read(), o_dofs1D, c_dofs1D);
|
||||||
@@ -6228,14 +6228,14 @@ static void PAHcurlVecH1IdentityApplyTranspose3D(const int c_dofs1D,
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PAHcurlVecH1IdentityApply2D(const int c_dofs1D,
|
void PAHcurlVecH1IdentityApply2D(const int c_dofs1D,
|
||||||
const int o_dofs1D,
|
const int o_dofs1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &Bclosed,
|
const Array<double> &Bclosed,
|
||||||
const Array<double> &Bopen,
|
const Array<double> &Bopen,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
auto Bc = Reshape(Bclosed.Read(), c_dofs1D, c_dofs1D);
|
auto Bc = Reshape(Bclosed.Read(), c_dofs1D, c_dofs1D);
|
||||||
auto Bo = Reshape(Bopen.Read(), o_dofs1D, c_dofs1D);
|
auto Bo = Reshape(Bopen.Read(), o_dofs1D, c_dofs1D);
|
||||||
@@ -6327,14 +6327,14 @@ static void PAHcurlVecH1IdentityApply2D(const int c_dofs1D,
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PAHcurlVecH1IdentityApplyTranspose2D(const int c_dofs1D,
|
void PAHcurlVecH1IdentityApplyTranspose2D(const int c_dofs1D,
|
||||||
const int o_dofs1D,
|
const int o_dofs1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &Bclosed,
|
const Array<double> &Bclosed,
|
||||||
const Array<double> &Bopen,
|
const Array<double> &Bopen,
|
||||||
const Vector &pa_data,
|
const Vector &pa_data,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
auto Bc = Reshape(Bclosed.Read(), c_dofs1D, c_dofs1D);
|
auto Bc = Reshape(Bclosed.Read(), c_dofs1D, c_dofs1D);
|
||||||
auto Bo = Reshape(Bopen.Read(), o_dofs1D, c_dofs1D);
|
auto Bo = Reshape(Bopen.Read(), o_dofs1D, c_dofs1D);
|
||||||
|
|||||||
+116
-116
@@ -539,12 +539,12 @@ void PAHdivMassApply3D(const int D1D,
|
|||||||
|
|
||||||
// PA H(div) div-div assemble 2D kernel
|
// PA H(div) div-div assemble 2D kernel
|
||||||
// NOTE: this is identical to PACurlCurlSetup3D
|
// NOTE: this is identical to PACurlCurlSetup3D
|
||||||
static void PADivDivSetup2D(const int Q1D,
|
void PADivDivSetup2D(const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
Vector &coeff_,
|
Vector &coeff_,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D;
|
const int NQ = Q1D*Q1D;
|
||||||
auto W = w.Read();
|
auto W = w.Read();
|
||||||
@@ -565,12 +565,12 @@ static void PADivDivSetup2D(const int Q1D,
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PADivDivSetup3D(const int Q1D,
|
void PADivDivSetup3D(const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
Vector &coeff_,
|
Vector &coeff_,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D*Q1D;
|
const int NQ = Q1D*Q1D*Q1D;
|
||||||
auto W = w.Read();
|
auto W = w.Read();
|
||||||
@@ -599,16 +599,16 @@ static void PADivDivSetup3D(const int Q1D,
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PADivDivApply2D(const int D1D,
|
void PADivDivApply2D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &Bo_,
|
const Array<double> &Bo_,
|
||||||
const Array<double> &Gc_,
|
const Array<double> &Gc_,
|
||||||
const Array<double> &Bot_,
|
const Array<double> &Bot_,
|
||||||
const Array<double> &Gct_,
|
const Array<double> &Gct_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
constexpr static int VDIM = 2;
|
constexpr static int VDIM = 2;
|
||||||
constexpr static int MAX_D1D = HDIV_MAX_D1D;
|
constexpr static int MAX_D1D = HDIV_MAX_D1D;
|
||||||
@@ -718,16 +718,16 @@ static void PADivDivApply2D(const int D1D,
|
|||||||
}); // end of element loop
|
}); // end of element loop
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PADivDivApply3D(const int D1D,
|
void PADivDivApply3D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &Bo_,
|
const Array<double> &Bo_,
|
||||||
const Array<double> &Gc_,
|
const Array<double> &Gc_,
|
||||||
const Array<double> &Bot_,
|
const Array<double> &Bot_,
|
||||||
const Array<double> &Gct_,
|
const Array<double> &Gct_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= HDIV_MAX_D1D, "Error: D1D > HDIV_MAX_D1D");
|
MFEM_VERIFY(D1D <= HDIV_MAX_D1D, "Error: D1D > HDIV_MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= HDIV_MAX_Q1D, "Error: Q1D > HDIV_MAX_Q1D");
|
MFEM_VERIFY(Q1D <= HDIV_MAX_Q1D, "Error: Q1D > HDIV_MAX_Q1D");
|
||||||
@@ -967,13 +967,13 @@ void DivDivIntegrator::AddMultPA(const Vector &x, Vector &y) const
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PADivDivAssembleDiagonal2D(const int D1D,
|
void PADivDivAssembleDiagonal2D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &Bo_,
|
const Array<double> &Bo_,
|
||||||
const Array<double> &Gc_,
|
const Array<double> &Gc_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
Vector &diag_)
|
Vector &diag_)
|
||||||
{
|
{
|
||||||
constexpr static int VDIM = 2;
|
constexpr static int VDIM = 2;
|
||||||
constexpr static int MAX_Q1D = HDIV_MAX_Q1D;
|
constexpr static int MAX_Q1D = HDIV_MAX_Q1D;
|
||||||
@@ -1023,13 +1023,13 @@ static void PADivDivAssembleDiagonal2D(const int D1D,
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PADivDivAssembleDiagonal3D(const int D1D,
|
void PADivDivAssembleDiagonal3D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &Bo_,
|
const Array<double> &Bo_,
|
||||||
const Array<double> &Gc_,
|
const Array<double> &Gc_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
Vector &diag_)
|
Vector &diag_)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= HDIV_MAX_D1D, "Error: D1D > HDIV_MAX_D1D");
|
MFEM_VERIFY(D1D <= HDIV_MAX_D1D, "Error: D1D > HDIV_MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= HDIV_MAX_Q1D, "Error: Q1D > HDIV_MAX_Q1D");
|
MFEM_VERIFY(Q1D <= HDIV_MAX_Q1D, "Error: Q1D > HDIV_MAX_Q1D");
|
||||||
@@ -1104,11 +1104,11 @@ void DivDivIntegrator::AssembleDiagonalPA(Vector& diag)
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA H(div)-L2 (div u, p) assemble 2D kernel
|
// PA H(div)-L2 (div u, p) assemble 2D kernel
|
||||||
static void PADivL2Setup2D(const int Q1D,
|
void PADivL2Setup2D(const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
Vector &coeff_,
|
Vector &coeff_,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D;
|
const int NQ = Q1D*Q1D;
|
||||||
auto W = w.Read();
|
auto W = w.Read();
|
||||||
@@ -1123,11 +1123,11 @@ static void PADivL2Setup2D(const int Q1D,
|
|||||||
});
|
});
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PADivL2Setup3D(const int Q1D,
|
void PADivL2Setup3D(const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
Vector &coeff_,
|
Vector &coeff_,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D*Q1D;
|
const int NQ = Q1D*Q1D*Q1D;
|
||||||
auto W = w.Read();
|
auto W = w.Read();
|
||||||
@@ -1225,16 +1225,16 @@ VectorFEDivergenceIntegrator::AssemblePA(const FiniteElementSpace &trial_fes,
|
|||||||
|
|
||||||
// Apply to x corresponding to DOF's in H(div) (trial), whose divergence is
|
// Apply to x corresponding to DOF's in H(div) (trial), whose divergence is
|
||||||
// integrated against L_2 test functions corresponding to y.
|
// integrated against L_2 test functions corresponding to y.
|
||||||
static void PAHdivL2Apply3D(const int D1D,
|
void PAHdivL2Apply3D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int L2D1D,
|
const int L2D1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &Bo_,
|
const Array<double> &Bo_,
|
||||||
const Array<double> &Gc_,
|
const Array<double> &Gc_,
|
||||||
const Array<double> &L2Bot_,
|
const Array<double> &L2Bot_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= HDIV_MAX_D1D, "Error: D1D > HDIV_MAX_D1D");
|
MFEM_VERIFY(D1D <= HDIV_MAX_D1D, "Error: D1D > HDIV_MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= HDIV_MAX_Q1D, "Error: Q1D > HDIV_MAX_Q1D");
|
MFEM_VERIFY(Q1D <= HDIV_MAX_Q1D, "Error: Q1D > HDIV_MAX_Q1D");
|
||||||
@@ -1388,16 +1388,16 @@ static void PAHdivL2Apply3D(const int D1D,
|
|||||||
|
|
||||||
// Apply to x corresponding to DOF's in H(div) (trial), whose divergence is
|
// Apply to x corresponding to DOF's in H(div) (trial), whose divergence is
|
||||||
// integrated against L_2 test functions corresponding to y.
|
// integrated against L_2 test functions corresponding to y.
|
||||||
static void PAHdivL2Apply2D(const int D1D,
|
void PAHdivL2Apply2D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int L2D1D,
|
const int L2D1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &Bo_,
|
const Array<double> &Bo_,
|
||||||
const Array<double> &Gc_,
|
const Array<double> &Gc_,
|
||||||
const Array<double> &L2Bot_,
|
const Array<double> &L2Bot_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
constexpr static int VDIM = 2;
|
constexpr static int VDIM = 2;
|
||||||
constexpr static int MAX_D1D = HDIV_MAX_D1D;
|
constexpr static int MAX_D1D = HDIV_MAX_D1D;
|
||||||
@@ -1494,16 +1494,16 @@ static void PAHdivL2Apply2D(const int D1D,
|
|||||||
}); // end of element loop
|
}); // end of element loop
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PAHdivL2ApplyTranspose3D(const int D1D,
|
void PAHdivL2ApplyTranspose3D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int L2D1D,
|
const int L2D1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &L2Bo_,
|
const Array<double> &L2Bo_,
|
||||||
const Array<double> &Gct_,
|
const Array<double> &Gct_,
|
||||||
const Array<double> &Bot_,
|
const Array<double> &Bot_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= HDIV_MAX_D1D, "Error: D1D > HDIV_MAX_D1D");
|
MFEM_VERIFY(D1D <= HDIV_MAX_D1D, "Error: D1D > HDIV_MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= HDIV_MAX_Q1D, "Error: Q1D > HDIV_MAX_Q1D");
|
MFEM_VERIFY(Q1D <= HDIV_MAX_Q1D, "Error: Q1D > HDIV_MAX_Q1D");
|
||||||
@@ -1656,16 +1656,16 @@ static void PAHdivL2ApplyTranspose3D(const int D1D,
|
|||||||
}); // end of element loop
|
}); // end of element loop
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PAHdivL2ApplyTranspose2D(const int D1D,
|
void PAHdivL2ApplyTranspose2D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int L2D1D,
|
const int L2D1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &L2Bo_,
|
const Array<double> &L2Bo_,
|
||||||
const Array<double> &Gct_,
|
const Array<double> &Gct_,
|
||||||
const Array<double> &Bot_,
|
const Array<double> &Bot_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_)
|
Vector &y_)
|
||||||
{
|
{
|
||||||
constexpr static int VDIM = 2;
|
constexpr static int VDIM = 2;
|
||||||
constexpr static int MAX_D1D = HDIV_MAX_D1D;
|
constexpr static int MAX_D1D = HDIV_MAX_D1D;
|
||||||
@@ -1791,16 +1791,16 @@ void VectorFEDivergenceIntegrator::AddMultTransposePA(const Vector &x,
|
|||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PAHdivL2AssembleDiagonal_ADAt_3D(const int D1D,
|
void PAHdivL2AssembleDiagonal_ADAt_3D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int L2D1D,
|
const int L2D1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &L2Bo_,
|
const Array<double> &L2Bo_,
|
||||||
const Array<double> &Gct_,
|
const Array<double> &Gct_,
|
||||||
const Array<double> &Bot_,
|
const Array<double> &Bot_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &D_,
|
const Vector &D_,
|
||||||
Vector &diag_)
|
Vector &diag_)
|
||||||
{
|
{
|
||||||
MFEM_VERIFY(D1D <= HDIV_MAX_D1D, "Error: D1D > HDIV_MAX_D1D");
|
MFEM_VERIFY(D1D <= HDIV_MAX_D1D, "Error: D1D > HDIV_MAX_D1D");
|
||||||
MFEM_VERIFY(Q1D <= HDIV_MAX_Q1D, "Error: Q1D > HDIV_MAX_Q1D");
|
MFEM_VERIFY(Q1D <= HDIV_MAX_Q1D, "Error: Q1D > HDIV_MAX_Q1D");
|
||||||
@@ -1916,16 +1916,16 @@ static void PAHdivL2AssembleDiagonal_ADAt_3D(const int D1D,
|
|||||||
}); // end of element loop
|
}); // end of element loop
|
||||||
}
|
}
|
||||||
|
|
||||||
static void PAHdivL2AssembleDiagonal_ADAt_2D(const int D1D,
|
void PAHdivL2AssembleDiagonal_ADAt_2D(const int D1D,
|
||||||
const int Q1D,
|
const int Q1D,
|
||||||
const int L2D1D,
|
const int L2D1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &L2Bo_,
|
const Array<double> &L2Bo_,
|
||||||
const Array<double> &Gct_,
|
const Array<double> &Gct_,
|
||||||
const Array<double> &Bot_,
|
const Array<double> &Bot_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &D_,
|
const Vector &D_,
|
||||||
Vector &diag_)
|
Vector &diag_)
|
||||||
{
|
{
|
||||||
constexpr static int VDIM = 2;
|
constexpr static int VDIM = 2;
|
||||||
|
|
||||||
|
|||||||
+21
-21
@@ -17,13 +17,13 @@ namespace mfem
|
|||||||
{
|
{
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EAMassAssemble1D(const int NE,
|
void EAMassAssemble1D(const int NE,
|
||||||
const Array<double> &basis,
|
const Array<double> &basis,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata,
|
Vector &eadata,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -67,13 +67,13 @@ static void EAMassAssemble1D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EAMassAssemble2D(const int NE,
|
void EAMassAssemble2D(const int NE,
|
||||||
const Array<double> &basis,
|
const Array<double> &basis,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata,
|
Vector &eadata,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -139,13 +139,13 @@ static void EAMassAssemble2D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void EAMassAssemble3D(const int NE,
|
void EAMassAssemble3D(const int NE,
|
||||||
const Array<double> &basis,
|
const Array<double> &basis,
|
||||||
const Vector &padata,
|
const Vector &padata,
|
||||||
Vector &eadata,
|
Vector &eadata,
|
||||||
const bool add,
|
const bool add,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
|
|||||||
+56
-56
@@ -155,12 +155,12 @@ void MassIntegrator::AssemblePA(const FiniteElementSpace &fes)
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PAMassAssembleDiagonal2D(const int NE,
|
void PAMassAssembleDiagonal2D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Vector &d,
|
const Vector &d,
|
||||||
Vector &y,
|
Vector &y,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -201,12 +201,12 @@ static void PAMassAssembleDiagonal2D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0>
|
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0>
|
||||||
static void SmemPAMassAssembleDiagonal2D(const int NE,
|
void SmemPAMassAssembleDiagonal2D(const int NE,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -267,12 +267,12 @@ static void SmemPAMassAssembleDiagonal2D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PAMassAssembleDiagonal3D(const int NE,
|
void PAMassAssembleDiagonal3D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Vector &d,
|
const Vector &d,
|
||||||
Vector &y,
|
Vector &y,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -336,12 +336,12 @@ static void PAMassAssembleDiagonal3D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void SmemPAMassAssembleDiagonal3D(const int NE,
|
void SmemPAMassAssembleDiagonal3D(const int NE,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -569,14 +569,14 @@ static void OccaPAMassApply3D(const int D1D,
|
|||||||
#endif // MFEM_USE_OCCA
|
#endif // MFEM_USE_OCCA
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PAMassApply2D(const int NE,
|
void PAMassApply2D(const int NE,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &bt_,
|
const Array<double> &bt_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -661,14 +661,14 @@ static void PAMassApply2D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0>
|
template<int T_D1D = 0, int T_Q1D = 0, int T_NBZ = 0>
|
||||||
static void SmemPAMassApply2D(const int NE,
|
void SmemPAMassApply2D(const int NE,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &bt_,
|
const Array<double> &bt_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
MFEM_CONTRACT_VAR(bt_);
|
MFEM_CONTRACT_VAR(bt_);
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
@@ -784,14 +784,14 @@ static void SmemPAMassApply2D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PAMassApply3D(const int NE,
|
void PAMassApply3D(const int NE,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &bt_,
|
const Array<double> &bt_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -925,14 +925,14 @@ static void PAMassApply3D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void SmemPAMassApply3D(const int NE,
|
void SmemPAMassApply3D(const int NE,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &bt_,
|
const Array<double> &bt_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
MFEM_CONTRACT_VAR(bt_);
|
MFEM_CONTRACT_VAR(bt_);
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
|
|||||||
@@ -22,12 +22,12 @@ namespace mfem
|
|||||||
// PA Vector Diffusion Integrator
|
// PA Vector Diffusion Integrator
|
||||||
|
|
||||||
// PA Diffusion Assemble 2D kernel
|
// PA Diffusion Assemble 2D kernel
|
||||||
static void PAVectorDiffusionSetup2D(const int Q1D,
|
void PAVectorDiffusionSetup2D(const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
const Vector &c,
|
const Vector &c,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D;
|
const int NQ = Q1D*Q1D;
|
||||||
auto W = w.Read();
|
auto W = w.Read();
|
||||||
@@ -59,12 +59,12 @@ static void PAVectorDiffusionSetup2D(const int Q1D,
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA Diffusion Assemble 3D kernel
|
// PA Diffusion Assemble 3D kernel
|
||||||
static void PAVectorDiffusionSetup3D(const int Q1D,
|
void PAVectorDiffusionSetup3D(const int Q1D,
|
||||||
const int NE,
|
const int NE,
|
||||||
const Array<double> &w,
|
const Array<double> &w,
|
||||||
const Vector &j,
|
const Vector &j,
|
||||||
const Vector &c,
|
const Vector &c,
|
||||||
Vector &op)
|
Vector &op)
|
||||||
{
|
{
|
||||||
const int NQ = Q1D*Q1D*Q1D;
|
const int NQ = Q1D*Q1D*Q1D;
|
||||||
auto W = w.Read();
|
auto W = w.Read();
|
||||||
@@ -251,7 +251,7 @@ void VectorDiffusionIntegrator::AssemblePA(const FiniteElementSpace &fes)
|
|||||||
}
|
}
|
||||||
|
|
||||||
// PA Diffusion Apply 2D kernel
|
// PA Diffusion Apply 2D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0, int T_VDIM = 0> static
|
template<int T_D1D = 0, int T_Q1D = 0, int T_VDIM = 0>
|
||||||
void PAVectorDiffusionApply2D(const int NE,
|
void PAVectorDiffusionApply2D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
@@ -374,7 +374,7 @@ void PAVectorDiffusionApply2D(const int NE,
|
|||||||
|
|
||||||
// PA Diffusion Apply 3D kernel
|
// PA Diffusion Apply 3D kernel
|
||||||
template<const int T_D1D = 0,
|
template<const int T_D1D = 0,
|
||||||
const int T_Q1D = 0> static
|
const int T_Q1D = 0>
|
||||||
void PAVectorDiffusionApply3D(const int NE,
|
void PAVectorDiffusionApply3D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
@@ -606,13 +606,13 @@ void VectorDiffusionIntegrator::AddMultPA(const Vector &x, Vector &y) const
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PAVectorDiffusionDiagonal2D(const int NE,
|
void PAVectorDiffusionDiagonal2D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Vector &d,
|
const Vector &d,
|
||||||
Vector &y,
|
Vector &y,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -673,13 +673,13 @@ static void PAVectorDiffusionDiagonal2D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PAVectorDiffusionDiagonal3D(const int NE,
|
void PAVectorDiffusionDiagonal3D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Vector &d,
|
const Vector &d,
|
||||||
Vector &y,
|
Vector &y,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
constexpr int DIM = 3;
|
constexpr int DIM = 3;
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
|
|||||||
+30
-30
@@ -104,14 +104,14 @@ void VectorMassIntegrator::AssemblePA(const FiniteElementSpace &fes)
|
|||||||
|
|
||||||
template<const int T_D1D = 0,
|
template<const int T_D1D = 0,
|
||||||
const int T_Q1D = 0>
|
const int T_Q1D = 0>
|
||||||
static void PAVectorMassApply2D(const int NE,
|
void PAVectorMassApply2D(const int NE,
|
||||||
const Array<double> &B_,
|
const Array<double> &B_,
|
||||||
const Array<double> &Bt_,
|
const Array<double> &Bt_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -201,14 +201,14 @@ static void PAVectorMassApply2D(const int NE,
|
|||||||
|
|
||||||
template<const int T_D1D = 0,
|
template<const int T_D1D = 0,
|
||||||
const int T_Q1D = 0>
|
const int T_Q1D = 0>
|
||||||
static void PAVectorMassApply3D(const int NE,
|
void PAVectorMassApply3D(const int NE,
|
||||||
const Array<double> &B_,
|
const Array<double> &B_,
|
||||||
const Array<double> &Bt_,
|
const Array<double> &Bt_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -379,13 +379,13 @@ void VectorMassIntegrator::AddMultPA(const Vector &x, Vector &y) const
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<const int T_D1D = 0, const int T_Q1D = 0>
|
template<const int T_D1D = 0, const int T_Q1D = 0>
|
||||||
static void PAVectorMassAssembleDiagonal2D(const int NE,
|
void PAVectorMassAssembleDiagonal2D(const int NE,
|
||||||
const Array<double> &B_,
|
const Array<double> &B_,
|
||||||
const Array<double> &Bt_,
|
const Array<double> &Bt_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
Vector &diag_,
|
Vector &diag_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -431,13 +431,13 @@ static void PAVectorMassAssembleDiagonal2D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<const int T_D1D = 0, const int T_Q1D = 0>
|
template<const int T_D1D = 0, const int T_Q1D = 0>
|
||||||
static void PAVectorMassAssembleDiagonal3D(const int NE,
|
void PAVectorMassAssembleDiagonal3D(const int NE,
|
||||||
const Array<double> &B_,
|
const Array<double> &B_,
|
||||||
const Array<double> &Bt_,
|
const Array<double> &Bt_,
|
||||||
const Vector &op_,
|
const Vector &op_,
|
||||||
Vector &diag_,
|
Vector &diag_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
|
|||||||
@@ -116,15 +116,15 @@ void VectorConvectionNLFIntegrator::AssemblePA(const FiniteElementSpace &fes)
|
|||||||
|
|
||||||
// PA Convection NL 2D kernel
|
// PA Convection NL 2D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PAConvectionNLApply2D(const int NE,
|
void PAConvectionNLApply2D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
const Vector &q_,
|
const Vector &q_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -252,15 +252,15 @@ static void PAConvectionNLApply2D(const int NE,
|
|||||||
|
|
||||||
// PA Convection NL 3D kernel
|
// PA Convection NL 3D kernel
|
||||||
template<int T_D1D = 0, int T_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0>
|
||||||
static void PAConvectionNLApply3D(const int NE,
|
void PAConvectionNLApply3D(const int NE,
|
||||||
const Array<double> &b,
|
const Array<double> &b,
|
||||||
const Array<double> &g,
|
const Array<double> &g,
|
||||||
const Array<double> &bt,
|
const Array<double> &bt,
|
||||||
const Vector &q_,
|
const Vector &q_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
constexpr int VDIM = 3;
|
constexpr int VDIM = 3;
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
@@ -558,14 +558,14 @@ static void PAConvectionNLApply3D(const int NE,
|
|||||||
}
|
}
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0, int T_MAX_D1D =0, int T_MAX_Q1D =0>
|
template<int T_D1D = 0, int T_Q1D = 0, int T_MAX_D1D =0, int T_MAX_Q1D =0>
|
||||||
static void SmemPAConvectionNLApply3D(const int NE,
|
void SmemPAConvectionNLApply3D(const int NE,
|
||||||
const Array<double> &b_,
|
const Array<double> &b_,
|
||||||
const Array<double> &g_,
|
const Array<double> &g_,
|
||||||
const Vector &d_,
|
const Vector &d_,
|
||||||
const Vector &x_,
|
const Vector &x_,
|
||||||
Vector &y_,
|
Vector &y_,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
constexpr int VDIM = 3;
|
constexpr int VDIM = 3;
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
|
|||||||
+17
-17
@@ -27,14 +27,14 @@ namespace quadrature_interpolator
|
|||||||
{
|
{
|
||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0, int MAX_D1D = 0, int MAX_Q1D = 0>
|
template<int T_D1D = 0, int T_Q1D = 0, int MAX_D1D = 0, int MAX_Q1D = 0>
|
||||||
static void Det2D(const int NE,
|
void Det2D(const int NE,
|
||||||
const double *b,
|
const double *b,
|
||||||
const double *g,
|
const double *g,
|
||||||
const double *x,
|
const double *x,
|
||||||
double *y,
|
double *y,
|
||||||
const int vdim = 1,
|
const int vdim = 1,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
constexpr int DIM = 2;
|
constexpr int DIM = 2;
|
||||||
static constexpr int NBZ = 1;
|
static constexpr int NBZ = 1;
|
||||||
@@ -79,15 +79,15 @@ static void Det2D(const int NE,
|
|||||||
|
|
||||||
template<int T_D1D = 0, int T_Q1D = 0, int MAX_D1D = 0, int MAX_Q1D = 0,
|
template<int T_D1D = 0, int T_Q1D = 0, int MAX_D1D = 0, int MAX_Q1D = 0,
|
||||||
bool SMEM = true>
|
bool SMEM = true>
|
||||||
static void Det3D(const int NE,
|
void Det3D(const int NE,
|
||||||
const double *b,
|
const double *b,
|
||||||
const double *g,
|
const double *g,
|
||||||
const double *x,
|
const double *x,
|
||||||
double *y,
|
double *y,
|
||||||
const int vdim = 1,
|
const int vdim = 1,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0,
|
const int q1d = 0,
|
||||||
Vector *d_buff = nullptr) // used only with SMEM = false
|
Vector *d_buff = nullptr) // used only with SMEM = false
|
||||||
{
|
{
|
||||||
constexpr int DIM = 3;
|
constexpr int DIM = 3;
|
||||||
static constexpr int MQ1 = T_Q1D ? T_Q1D : MAX_Q1D;
|
static constexpr int MQ1 = T_Q1D ? T_Q1D : MAX_Q1D;
|
||||||
|
|||||||
+14
-14
@@ -31,13 +31,13 @@ namespace quadrature_interpolator
|
|||||||
template<QVectorLayout Q_LAYOUT,
|
template<QVectorLayout Q_LAYOUT,
|
||||||
int T_VDIM = 0, int T_D1D = 0, int T_Q1D = 0,
|
int T_VDIM = 0, int T_D1D = 0, int T_Q1D = 0,
|
||||||
int T_NBZ = 1, int MAX_D1D = 0, int MAX_Q1D = 0>
|
int T_NBZ = 1, int MAX_D1D = 0, int MAX_Q1D = 0>
|
||||||
static void Values2D(const int NE,
|
void Values2D(const int NE,
|
||||||
const double *b_,
|
const double *b_,
|
||||||
const double *x_,
|
const double *x_,
|
||||||
double *y_,
|
double *y_,
|
||||||
const int vdim = 0,
|
const int vdim = 0,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
static constexpr int NBZ = T_NBZ ? T_NBZ : 1;
|
static constexpr int NBZ = T_NBZ ? T_NBZ : 1;
|
||||||
|
|
||||||
@@ -95,13 +95,13 @@ static void Values2D(const int NE,
|
|||||||
template<QVectorLayout Q_LAYOUT,
|
template<QVectorLayout Q_LAYOUT,
|
||||||
int T_VDIM = 0, int T_D1D = 0, int T_Q1D = 0,
|
int T_VDIM = 0, int T_D1D = 0, int T_Q1D = 0,
|
||||||
int MAX_D1D = 0, int MAX_Q1D = 0>
|
int MAX_D1D = 0, int MAX_Q1D = 0>
|
||||||
static void Values3D(const int NE,
|
void Values3D(const int NE,
|
||||||
const double *b_,
|
const double *b_,
|
||||||
const double *x_,
|
const double *x_,
|
||||||
double *y_,
|
double *y_,
|
||||||
const int vdim = 0,
|
const int vdim = 0,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
|
|||||||
+18
-18
@@ -31,15 +31,15 @@ namespace quadrature_interpolator
|
|||||||
template<QVectorLayout Q_LAYOUT, bool GRAD_PHYS,
|
template<QVectorLayout Q_LAYOUT, bool GRAD_PHYS,
|
||||||
int T_VDIM = 0, int T_D1D = 0, int T_Q1D = 0,
|
int T_VDIM = 0, int T_D1D = 0, int T_Q1D = 0,
|
||||||
int T_NBZ = 1, int MAX_D1D = 0, int MAX_Q1D = 0>
|
int T_NBZ = 1, int MAX_D1D = 0, int MAX_Q1D = 0>
|
||||||
static void Derivatives2D(const int NE,
|
void Derivatives2D(const int NE,
|
||||||
const double *b_,
|
const double *b_,
|
||||||
const double *g_,
|
const double *g_,
|
||||||
const double *j_,
|
const double *j_,
|
||||||
const double *x_,
|
const double *x_,
|
||||||
double *y_,
|
double *y_,
|
||||||
const int vdim = 0,
|
const int vdim = 0,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
@@ -139,15 +139,15 @@ static void Derivatives2D(const int NE,
|
|||||||
template<QVectorLayout Q_LAYOUT, bool GRAD_PHYS,
|
template<QVectorLayout Q_LAYOUT, bool GRAD_PHYS,
|
||||||
int T_VDIM = 0, int T_D1D = 0, int T_Q1D = 0,
|
int T_VDIM = 0, int T_D1D = 0, int T_Q1D = 0,
|
||||||
int MAX_D1D = 0, int MAX_Q1D = 0>
|
int MAX_D1D = 0, int MAX_Q1D = 0>
|
||||||
static void Derivatives3D(const int NE,
|
void Derivatives3D(const int NE,
|
||||||
const double *b_,
|
const double *b_,
|
||||||
const double *g_,
|
const double *g_,
|
||||||
const double *j_,
|
const double *j_,
|
||||||
const double *x_,
|
const double *x_,
|
||||||
double *y_,
|
double *y_,
|
||||||
const int vdim = 0,
|
const int vdim = 0,
|
||||||
const int d1d = 0,
|
const int d1d = 0,
|
||||||
const int q1d = 0)
|
const int q1d = 0)
|
||||||
{
|
{
|
||||||
const int D1D = T_D1D ? T_D1D : d1d;
|
const int D1D = T_D1D ? T_D1D : d1d;
|
||||||
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
const int Q1D = T_Q1D ? T_Q1D : q1d;
|
||||||
|
|||||||
+20
-20
@@ -61,16 +61,16 @@ namespace quadrature_interpolator
|
|||||||
// * assumes 'e_vec' is using ElementDofOrdering::NATIVE,
|
// * assumes 'e_vec' is using ElementDofOrdering::NATIVE,
|
||||||
// * assumes 'maps.mode == FULL'.
|
// * assumes 'maps.mode == FULL'.
|
||||||
template<const int T_VDIM, const int T_ND, const int T_NQ>
|
template<const int T_VDIM, const int T_ND, const int T_NQ>
|
||||||
static void Eval2D(const int NE,
|
void Eval2D(const int NE,
|
||||||
const int vdim,
|
const int vdim,
|
||||||
const QVectorLayout q_layout,
|
const QVectorLayout q_layout,
|
||||||
const GeometricFactors *geom,
|
const GeometricFactors *geom,
|
||||||
const DofToQuad &maps,
|
const DofToQuad &maps,
|
||||||
const Vector &e_vec,
|
const Vector &e_vec,
|
||||||
Vector &q_val,
|
Vector &q_val,
|
||||||
Vector &q_der,
|
Vector &q_der,
|
||||||
Vector &q_det,
|
Vector &q_det,
|
||||||
const int eval_flags)
|
const int eval_flags)
|
||||||
{
|
{
|
||||||
using QI = QuadratureInterpolator;
|
using QI = QuadratureInterpolator;
|
||||||
|
|
||||||
@@ -209,16 +209,16 @@ static void Eval2D(const int NE,
|
|||||||
// * assumes 'e_vec' is using ElementDofOrdering::NATIVE,
|
// * assumes 'e_vec' is using ElementDofOrdering::NATIVE,
|
||||||
// * assumes 'maps.mode == FULL'.
|
// * assumes 'maps.mode == FULL'.
|
||||||
template<const int T_VDIM, const int T_ND, const int T_NQ>
|
template<const int T_VDIM, const int T_ND, const int T_NQ>
|
||||||
static void Eval3D(const int NE,
|
void Eval3D(const int NE,
|
||||||
const int vdim,
|
const int vdim,
|
||||||
const QVectorLayout q_layout,
|
const QVectorLayout q_layout,
|
||||||
const GeometricFactors *geom,
|
const GeometricFactors *geom,
|
||||||
const DofToQuad &maps,
|
const DofToQuad &maps,
|
||||||
const Vector &e_vec,
|
const Vector &e_vec,
|
||||||
Vector &q_val,
|
Vector &q_val,
|
||||||
Vector &q_der,
|
Vector &q_der,
|
||||||
Vector &q_det,
|
Vector &q_det,
|
||||||
const int eval_flags)
|
const int eval_flags)
|
||||||
{
|
{
|
||||||
using QI = QuadratureInterpolator;
|
using QI = QuadratureInterpolator;
|
||||||
|
|
||||||
|
|||||||
+1
-1
@@ -1306,7 +1306,7 @@ namespace internal
|
|||||||
// MFEM_FORALL-based copy kernel -- used by protected methods below.
|
// MFEM_FORALL-based copy kernel -- used by protected methods below.
|
||||||
// Needed as a workaround for the nvcc restriction that methods with MFEM_FORALL
|
// Needed as a workaround for the nvcc restriction that methods with MFEM_FORALL
|
||||||
// in them must to be public.
|
// in them must to be public.
|
||||||
static inline void device_copy(double *d_dest, const double *d_src, int size)
|
inline void device_copy(double *d_dest, const double *d_src, int size)
|
||||||
{
|
{
|
||||||
MFEM_FORALL(i, size, d_dest[i] = d_src[i];);
|
MFEM_FORALL(i, size, d_dest[i] = d_src[i];);
|
||||||
}
|
}
|
||||||
|
|||||||
+1
-1
@@ -234,7 +234,7 @@ void TMOPRefinerEstimator::SetTriIntRules()
|
|||||||
|
|
||||||
// Reftype = 0 // original element
|
// Reftype = 0 // original element
|
||||||
const int Nvert = 3, NEsplit = 1;
|
const int Nvert = 3, NEsplit = 1;
|
||||||
Mesh meshsplit(2, Nvert, NEsplit, 0 ,2);
|
Mesh meshsplit(2, Nvert, NEsplit, 0,2);
|
||||||
const double tri_v[3][2] =
|
const double tri_v[3][2] =
|
||||||
{
|
{
|
||||||
{0, 0}, {1, 0}, {0, 1}
|
{0, 0}, {1, 0}, {0, 1}
|
||||||
|
|||||||
+7
-1
@@ -164,7 +164,13 @@ __device__ void abort_msg(T & msg)
|
|||||||
#endif
|
#endif
|
||||||
|
|
||||||
// Abort inside a device kernel
|
// Abort inside a device kernel
|
||||||
#if defined(__CUDA_ARCH__)
|
#if defined(__CUDA_ARCH__) && defined(_WIN32)
|
||||||
|
#define MFEM_ABORT_KERNEL(msg) \
|
||||||
|
{ \
|
||||||
|
printf(msg); \
|
||||||
|
__debugbreak(); \
|
||||||
|
}
|
||||||
|
#elif defined(__CUDA_ARCH__)
|
||||||
#define MFEM_ABORT_KERNEL(msg) \
|
#define MFEM_ABORT_KERNEL(msg) \
|
||||||
{ \
|
{ \
|
||||||
printf(msg); \
|
printf(msg); \
|
||||||
|
|||||||
@@ -79,7 +79,13 @@ int isockstream::establish()
|
|||||||
int on=1;
|
int on=1;
|
||||||
setsockopt(port, SOL_SOCKET, SO_REUSEADDR, (char *)(&on), sizeof(on));
|
setsockopt(port, SOL_SOCKET, SO_REUSEADDR, (char *)(&on), sizeof(on));
|
||||||
|
|
||||||
if (bind(port,(const sockaddr*)&sa,(socklen_t)sizeof(struct sockaddr_in)) < 0)
|
if (bind(
|
||||||
|
#ifdef _WIN32
|
||||||
|
(SOCKET)port
|
||||||
|
#else
|
||||||
|
port
|
||||||
|
#endif
|
||||||
|
,(const sockaddr*)&sa,(socklen_t)sizeof(struct sockaddr_in)) < 0)
|
||||||
{
|
{
|
||||||
mfem::err << "isockstream::establish(): bind() failed!" << endl;
|
mfem::err << "isockstream::establish(): bind() failed!" << endl;
|
||||||
close(port);
|
close(port);
|
||||||
|
|||||||
+3325
-2996
File diff suppressed because it is too large
Load Diff
+1
-1
@@ -378,7 +378,7 @@ void Swap(T &a, T &b)
|
|||||||
b = tmp;
|
b = tmp;
|
||||||
}
|
}
|
||||||
|
|
||||||
const double Epsilon = std::numeric_limits<double>::epsilon();
|
constexpr double Epsilon = std::numeric_limits<double>::epsilon();
|
||||||
|
|
||||||
/// Utility function used in CalcSingularvalue<3>.
|
/// Utility function used in CalcSingularvalue<3>.
|
||||||
MFEM_HOST_DEVICE static inline
|
MFEM_HOST_DEVICE static inline
|
||||||
|
|||||||
+1
-1
@@ -1677,7 +1677,7 @@ int KINSolver::GradientMult(N_Vector v, N_Vector Jv, N_Vector u,
|
|||||||
|
|
||||||
// Wrapper for evaluating linear systems J u = b
|
// Wrapper for evaluating linear systems J u = b
|
||||||
int KINSolver::LinSysSetup(N_Vector u, N_Vector, SUNMatrix J,
|
int KINSolver::LinSysSetup(N_Vector u, N_Vector, SUNMatrix J,
|
||||||
void *, N_Vector , N_Vector )
|
void *, N_Vector, N_Vector )
|
||||||
{
|
{
|
||||||
const SundialsNVector mfem_u(u);
|
const SundialsNVector mfem_u(u);
|
||||||
KINSolver *self = static_cast<KINSolver*>(GET_CONTENT(J));
|
KINSolver *self = static_cast<KINSolver*>(GET_CONTENT(J));
|
||||||
|
|||||||
+27
-27
@@ -2850,14 +2850,14 @@ void Mesh::Make3D(int nx, int ny, int nz, Element::Type type,
|
|||||||
y = sfc[3*k + 1];
|
y = sfc[3*k + 1];
|
||||||
z = sfc[3*k + 2];
|
z = sfc[3*k + 2];
|
||||||
|
|
||||||
ind[0] = VTX(x , y , z );
|
ind[0] = VTX(x, y, z );
|
||||||
ind[1] = VTX(x+1, y , z );
|
ind[1] = VTX(x+1, y, z );
|
||||||
ind[2] = VTX(x+1, y+1, z );
|
ind[2] = VTX(x+1, y+1, z );
|
||||||
ind[3] = VTX(x , y+1, z );
|
ind[3] = VTX(x, y+1, z );
|
||||||
ind[4] = VTX(x , y , z+1);
|
ind[4] = VTX(x, y, z+1);
|
||||||
ind[5] = VTX(x+1, y , z+1);
|
ind[5] = VTX(x+1, y, z+1);
|
||||||
ind[6] = VTX(x+1, y+1, z+1);
|
ind[6] = VTX(x+1, y+1, z+1);
|
||||||
ind[7] = VTX(x , y+1, z+1);
|
ind[7] = VTX(x, y+1, z+1);
|
||||||
|
|
||||||
AddHex(ind, 1);
|
AddHex(ind, 1);
|
||||||
}
|
}
|
||||||
@@ -2870,12 +2870,12 @@ void Mesh::Make3D(int nx, int ny, int nz, Element::Type type,
|
|||||||
{
|
{
|
||||||
for (x = 0; x < nx; x++)
|
for (x = 0; x < nx; x++)
|
||||||
{
|
{
|
||||||
ind[0] = VTX(x , y , z );
|
ind[0] = VTX(x, y, z );
|
||||||
ind[1] = VTX(x+1, y , z );
|
ind[1] = VTX(x+1, y, z );
|
||||||
ind[2] = VTX(x+1, y+1, z );
|
ind[2] = VTX(x+1, y+1, z );
|
||||||
ind[3] = VTX(x , y+1, z );
|
ind[3] = VTX(x, y+1, z );
|
||||||
ind[4] = VTX(x , y , z+1);
|
ind[4] = VTX(x, y, z+1);
|
||||||
ind[5] = VTX(x+1, y , z+1);
|
ind[5] = VTX(x+1, y, z+1);
|
||||||
ind[6] = VTX(x+1, y+1, z+1);
|
ind[6] = VTX(x+1, y+1, z+1);
|
||||||
ind[7] = VTX( x, y+1, z+1);
|
ind[7] = VTX( x, y+1, z+1);
|
||||||
if (type == Element::TETRAHEDRON)
|
if (type == Element::TETRAHEDRON)
|
||||||
@@ -2906,10 +2906,10 @@ void Mesh::Make3D(int nx, int ny, int nz, Element::Type type,
|
|||||||
{
|
{
|
||||||
for (x = 0; x < nx; x++)
|
for (x = 0; x < nx; x++)
|
||||||
{
|
{
|
||||||
ind[0] = VTX(x , y , 0);
|
ind[0] = VTX(x, y, 0);
|
||||||
ind[1] = VTX(x , y+1, 0);
|
ind[1] = VTX(x, y+1, 0);
|
||||||
ind[2] = VTX(x+1, y+1, 0);
|
ind[2] = VTX(x+1, y+1, 0);
|
||||||
ind[3] = VTX(x+1, y , 0);
|
ind[3] = VTX(x+1, y, 0);
|
||||||
if (type == Element::TETRAHEDRON)
|
if (type == Element::TETRAHEDRON)
|
||||||
{
|
{
|
||||||
AddBdrQuadAsTriangles(ind, 1);
|
AddBdrQuadAsTriangles(ind, 1);
|
||||||
@@ -2929,10 +2929,10 @@ void Mesh::Make3D(int nx, int ny, int nz, Element::Type type,
|
|||||||
{
|
{
|
||||||
for (x = 0; x < nx; x++)
|
for (x = 0; x < nx; x++)
|
||||||
{
|
{
|
||||||
ind[0] = VTX(x , y , nz);
|
ind[0] = VTX(x, y, nz);
|
||||||
ind[1] = VTX(x+1, y , nz);
|
ind[1] = VTX(x+1, y, nz);
|
||||||
ind[2] = VTX(x+1, y+1, nz);
|
ind[2] = VTX(x+1, y+1, nz);
|
||||||
ind[3] = VTX(x , y+1, nz);
|
ind[3] = VTX(x, y+1, nz);
|
||||||
if (type == Element::TETRAHEDRON)
|
if (type == Element::TETRAHEDRON)
|
||||||
{
|
{
|
||||||
AddBdrQuadAsTriangles(ind, 6);
|
AddBdrQuadAsTriangles(ind, 6);
|
||||||
@@ -2952,10 +2952,10 @@ void Mesh::Make3D(int nx, int ny, int nz, Element::Type type,
|
|||||||
{
|
{
|
||||||
for (y = 0; y < ny; y++)
|
for (y = 0; y < ny; y++)
|
||||||
{
|
{
|
||||||
ind[0] = VTX(0 , y , z );
|
ind[0] = VTX(0, y, z );
|
||||||
ind[1] = VTX(0 , y , z+1);
|
ind[1] = VTX(0, y, z+1);
|
||||||
ind[2] = VTX(0 , y+1, z+1);
|
ind[2] = VTX(0, y+1, z+1);
|
||||||
ind[3] = VTX(0 , y+1, z );
|
ind[3] = VTX(0, y+1, z );
|
||||||
if (type == Element::TETRAHEDRON)
|
if (type == Element::TETRAHEDRON)
|
||||||
{
|
{
|
||||||
AddBdrQuadAsTriangles(ind, 5);
|
AddBdrQuadAsTriangles(ind, 5);
|
||||||
@@ -2971,10 +2971,10 @@ void Mesh::Make3D(int nx, int ny, int nz, Element::Type type,
|
|||||||
{
|
{
|
||||||
for (y = 0; y < ny; y++)
|
for (y = 0; y < ny; y++)
|
||||||
{
|
{
|
||||||
ind[0] = VTX(nx, y , z );
|
ind[0] = VTX(nx, y, z );
|
||||||
ind[1] = VTX(nx, y+1, z );
|
ind[1] = VTX(nx, y+1, z );
|
||||||
ind[2] = VTX(nx, y+1, z+1);
|
ind[2] = VTX(nx, y+1, z+1);
|
||||||
ind[3] = VTX(nx, y , z+1);
|
ind[3] = VTX(nx, y, z+1);
|
||||||
if (type == Element::TETRAHEDRON)
|
if (type == Element::TETRAHEDRON)
|
||||||
{
|
{
|
||||||
AddBdrQuadAsTriangles(ind, 3);
|
AddBdrQuadAsTriangles(ind, 3);
|
||||||
@@ -2990,10 +2990,10 @@ void Mesh::Make3D(int nx, int ny, int nz, Element::Type type,
|
|||||||
{
|
{
|
||||||
for (z = 0; z < nz; z++)
|
for (z = 0; z < nz; z++)
|
||||||
{
|
{
|
||||||
ind[0] = VTX(x , 0, z );
|
ind[0] = VTX(x, 0, z );
|
||||||
ind[1] = VTX(x+1, 0, z );
|
ind[1] = VTX(x+1, 0, z );
|
||||||
ind[2] = VTX(x+1, 0, z+1);
|
ind[2] = VTX(x+1, 0, z+1);
|
||||||
ind[3] = VTX(x , 0, z+1);
|
ind[3] = VTX(x, 0, z+1);
|
||||||
if (type == Element::TETRAHEDRON)
|
if (type == Element::TETRAHEDRON)
|
||||||
{
|
{
|
||||||
AddBdrQuadAsTriangles(ind, 2);
|
AddBdrQuadAsTriangles(ind, 2);
|
||||||
@@ -3009,8 +3009,8 @@ void Mesh::Make3D(int nx, int ny, int nz, Element::Type type,
|
|||||||
{
|
{
|
||||||
for (z = 0; z < nz; z++)
|
for (z = 0; z < nz; z++)
|
||||||
{
|
{
|
||||||
ind[0] = VTX(x , ny, z );
|
ind[0] = VTX(x, ny, z );
|
||||||
ind[1] = VTX(x , ny, z+1);
|
ind[1] = VTX(x, ny, z+1);
|
||||||
ind[2] = VTX(x+1, ny, z+1);
|
ind[2] = VTX(x+1, ny, z+1);
|
||||||
ind[3] = VTX(x+1, ny, z );
|
ind[3] = VTX(x+1, ny, z );
|
||||||
if (type == Element::TETRAHEDRON)
|
if (type == Element::TETRAHEDRON)
|
||||||
|
|||||||
+1
-1
@@ -612,7 +612,7 @@ int NCMesh::NewSegment(int n0, int n1, int attr, int vattr1, int vattr2)
|
|||||||
// get (degenerate) faces and assign face attributes
|
// get (degenerate) faces and assign face attributes
|
||||||
int v0 = el.node[0], v1 = el.node[1];
|
int v0 = el.node[0], v1 = el.node[1];
|
||||||
faces.Get(v0, v0, v0, v0)->attribute = vattr1;
|
faces.Get(v0, v0, v0, v0)->attribute = vattr1;
|
||||||
faces.Get(v1, v1, v1 ,v1)->attribute = vattr2;
|
faces.Get(v1, v1, v1,v1)->attribute = vattr2;
|
||||||
|
|
||||||
return new_id;
|
return new_id;
|
||||||
}
|
}
|
||||||
|
|||||||
@@ -1159,8 +1159,8 @@ static double u0(const Vector &x) { return sin(3.0 * PI * (x[1] + x[0])); }
|
|||||||
|
|
||||||
enum {NORM, AREA};
|
enum {NORM, AREA};
|
||||||
|
|
||||||
static double qf(const int order, const int ker, Mesh &m,
|
double qf(const int order, const int ker, Mesh &m,
|
||||||
FiniteElementSpace &fes, GridFunction &u)
|
FiniteElementSpace &fes, GridFunction &u)
|
||||||
{
|
{
|
||||||
const Geometry::Type type = m.GetElementBaseGeometry(0);
|
const Geometry::Type type = m.GetElementBaseGeometry(0);
|
||||||
const IntegrationRule &ir(IntRules.Get(type, order));
|
const IntegrationRule &ir(IntRules.Get(type, order));
|
||||||
|
|||||||
@@ -483,11 +483,11 @@ void ScreenedPoisson::AssembleElementVector(const FiniteElement &el,
|
|||||||
pval=shapef*elfun;
|
pval=shapef*elfun;
|
||||||
if (fval>0.0)
|
if (fval>0.0)
|
||||||
{
|
{
|
||||||
elvect.Add( -w , shapef);
|
elvect.Add( -w, shapef);
|
||||||
}
|
}
|
||||||
else if (fval<0.0)
|
else if (fval<0.0)
|
||||||
{
|
{
|
||||||
elvect.Add( w , shapef);
|
elvect.Add( w, shapef);
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
@@ -523,7 +523,7 @@ void ScreenedPoisson::AssembleElementGrad(const FiniteElement &el,
|
|||||||
el.CalcPhysDShape(trans, B);
|
el.CalcPhysDShape(trans, B);
|
||||||
el.CalcPhysShape(trans,shapef);
|
el.CalcPhysShape(trans,shapef);
|
||||||
|
|
||||||
AddMult_a_VVt(w , shapef, elmat);
|
AddMult_a_VVt(w, shapef, elmat);
|
||||||
AddMult_a_AAt(w * diffcoef, B, elmat);
|
AddMult_a_AAt(w * diffcoef, B, elmat);
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
@@ -683,11 +683,11 @@ void PUMPLaplacian::AssembleElementVector(const FiniteElement &el,
|
|||||||
// add the external load -1 if tval > 0.0; 1 if tval < 0.0;
|
// add the external load -1 if tval > 0.0; 1 if tval < 0.0;
|
||||||
if (tval>0.0)
|
if (tval>0.0)
|
||||||
{
|
{
|
||||||
elvect.Add( -w*fval , shapef);
|
elvect.Add( -w*fval, shapef);
|
||||||
}
|
}
|
||||||
else if (tval<0.0)
|
else if (tval<0.0)
|
||||||
{
|
{
|
||||||
elvect.Add( w*fval , shapef);
|
elvect.Add( w*fval, shapef);
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|||||||
@@ -117,7 +117,7 @@ int main(int argc, char *argv[])
|
|||||||
MPI_Finalize();
|
MPI_Finalize();
|
||||||
return 1;
|
return 1;
|
||||||
}
|
}
|
||||||
if (prob >3 || prob <0) prob = 0; // default problem = H1
|
if (prob >3 || prob <0) { prob = 0; } // default problem = H1
|
||||||
if (prob == 3)
|
if (prob == 3)
|
||||||
{
|
{
|
||||||
if (kappa < 0)
|
if (kappa < 0)
|
||||||
@@ -221,7 +221,7 @@ int main(int argc, char *argv[])
|
|||||||
gradu = new VectorFunctionCoefficient(dim,gradu_exact);
|
gradu = new VectorFunctionCoefficient(dim,gradu_exact);
|
||||||
b.AddDomainIntegrator(new DomainLFIntegrator(*f));
|
b.AddDomainIntegrator(new DomainLFIntegrator(*f));
|
||||||
b.AddBdrFaceIntegrator(
|
b.AddBdrFaceIntegrator(
|
||||||
new DGDirichletLFIntegrator(*scalar_u, one, sigma, kappa));
|
new DGDirichletLFIntegrator(*scalar_u, one, sigma, kappa));
|
||||||
a.AddDomainIntegrator(new DiffusionIntegrator(one));
|
a.AddDomainIntegrator(new DiffusionIntegrator(one));
|
||||||
a.AddInteriorFaceIntegrator(new DGDiffusionIntegrator(one, sigma, kappa));
|
a.AddInteriorFaceIntegrator(new DGDiffusionIntegrator(one, sigma, kappa));
|
||||||
a.AddBdrFaceIntegrator(new DGDiffusionIntegrator(one, sigma, kappa));
|
a.AddBdrFaceIntegrator(new DGDiffusionIntegrator(one, sigma, kappa));
|
||||||
@@ -291,8 +291,8 @@ int main(int argc, char *argv[])
|
|||||||
|
|
||||||
x = *X;
|
x = *X;
|
||||||
JumpScaling js(1.0, jump_scaling_type == 2 ? JumpScaling::P_SQUARED_OVER_H
|
JumpScaling js(1.0, jump_scaling_type == 2 ? JumpScaling::P_SQUARED_OVER_H
|
||||||
: jump_scaling_type == 1 ? JumpScaling::ONE_OVER_H
|
: jump_scaling_type == 1 ? JumpScaling::ONE_OVER_H
|
||||||
: JumpScaling::CONSTANT);
|
: JumpScaling::CONSTANT);
|
||||||
switch (prob)
|
switch (prob)
|
||||||
{
|
{
|
||||||
case 0: rates.AddH1GridFunction(&x,scalar_u,gradu); break;
|
case 0: rates.AddH1GridFunction(&x,scalar_u,gradu); break;
|
||||||
@@ -305,7 +305,7 @@ int main(int argc, char *argv[])
|
|||||||
delete B;
|
delete B;
|
||||||
delete A;
|
delete A;
|
||||||
|
|
||||||
if (l==pr) break;
|
if (l==pr) { break; }
|
||||||
|
|
||||||
pmesh->UniformRefinement();
|
pmesh->UniformRefinement();
|
||||||
fespace->Update();
|
fespace->Update();
|
||||||
|
|||||||
@@ -104,7 +104,7 @@ int main(int argc, char *argv[])
|
|||||||
args.PrintUsage(cout);
|
args.PrintUsage(cout);
|
||||||
return 1;
|
return 1;
|
||||||
}
|
}
|
||||||
if (prob >3 || prob <0) prob = 0; // default problem = H1
|
if (prob >3 || prob <0) { prob = 0; } // default problem = H1
|
||||||
if (prob == 3)
|
if (prob == 3)
|
||||||
{
|
{
|
||||||
if (kappa < 0)
|
if (kappa < 0)
|
||||||
@@ -194,7 +194,7 @@ int main(int argc, char *argv[])
|
|||||||
gradu = new VectorFunctionCoefficient(dim,gradu_exact);
|
gradu = new VectorFunctionCoefficient(dim,gradu_exact);
|
||||||
b.AddDomainIntegrator(new DomainLFIntegrator(*f));
|
b.AddDomainIntegrator(new DomainLFIntegrator(*f));
|
||||||
b.AddBdrFaceIntegrator(
|
b.AddBdrFaceIntegrator(
|
||||||
new DGDirichletLFIntegrator(*scalar_u, one, sigma, kappa));
|
new DGDirichletLFIntegrator(*scalar_u, one, sigma, kappa));
|
||||||
a.AddDomainIntegrator(new DiffusionIntegrator(one));
|
a.AddDomainIntegrator(new DiffusionIntegrator(one));
|
||||||
a.AddInteriorFaceIntegrator(new DGDiffusionIntegrator(one, sigma, kappa));
|
a.AddInteriorFaceIntegrator(new DGDiffusionIntegrator(one, sigma, kappa));
|
||||||
a.AddBdrFaceIntegrator(new DGDiffusionIntegrator(one, sigma, kappa));
|
a.AddBdrFaceIntegrator(new DGDiffusionIntegrator(one, sigma, kappa));
|
||||||
@@ -224,8 +224,8 @@ int main(int argc, char *argv[])
|
|||||||
}
|
}
|
||||||
|
|
||||||
JumpScaling js(1.0, jump_scaling_type == 2 ? JumpScaling::P_SQUARED_OVER_H
|
JumpScaling js(1.0, jump_scaling_type == 2 ? JumpScaling::P_SQUARED_OVER_H
|
||||||
: jump_scaling_type == 1 ? JumpScaling::ONE_OVER_H
|
: jump_scaling_type == 1 ? JumpScaling::ONE_OVER_H
|
||||||
: JumpScaling::CONSTANT);
|
: JumpScaling::CONSTANT);
|
||||||
|
|
||||||
switch (prob)
|
switch (prob)
|
||||||
{
|
{
|
||||||
@@ -235,7 +235,7 @@ int main(int argc, char *argv[])
|
|||||||
case 3: rates.AddL2GridFunction(&x,scalar_u,gradu,&one,js); break;
|
case 3: rates.AddL2GridFunction(&x,scalar_u,gradu,&one,js); break;
|
||||||
}
|
}
|
||||||
|
|
||||||
if (l==sr) break;
|
if (l==sr) { break; }
|
||||||
|
|
||||||
mesh->UniformRefinement();
|
mesh->UniformRefinement();
|
||||||
fespace->Update();
|
fespace->Update();
|
||||||
|
|||||||
@@ -98,7 +98,7 @@ TEST_CASE("DoF Transformation Classes",
|
|||||||
double uAv = A.InnerProduct(v, u);
|
double uAv = A.InnerProduct(v, u);
|
||||||
|
|
||||||
REQUIRE(fabs(uAv - At.InnerProduct(vt, u )) < tol * fabs(uAv));
|
REQUIRE(fabs(uAv - At.InnerProduct(vt, u )) < tol * fabs(uAv));
|
||||||
REQUIRE(fabs(uAv - tA.InnerProduct(v , ut)) < tol * fabs(uAv));
|
REQUIRE(fabs(uAv - tA.InnerProduct(v, ut)) < tol * fabs(uAv));
|
||||||
REQUIRE(fabs(uAv - tAt.InnerProduct(vt, ut)) < tol * fabs(uAv));
|
REQUIRE(fabs(uAv - tAt.InnerProduct(vt, ut)) < tol * fabs(uAv));
|
||||||
}
|
}
|
||||||
SECTION("Inner product of a primal vector and a dual vector")
|
SECTION("Inner product of a primal vector and a dual vector")
|
||||||
@@ -119,7 +119,7 @@ TEST_CASE("DoF Transformation Classes",
|
|||||||
double fAv = A.InnerProduct(v, f);
|
double fAv = A.InnerProduct(v, f);
|
||||||
|
|
||||||
REQUIRE(fabs(fAv - At.InnerProduct(vt, f )) < tol * fabs(fAv));
|
REQUIRE(fabs(fAv - At.InnerProduct(vt, f )) < tol * fabs(fAv));
|
||||||
REQUIRE(fabs(fAv - tA.InnerProduct(v , ft)) < tol * fabs(fAv));
|
REQUIRE(fabs(fAv - tA.InnerProduct(v, ft)) < tol * fabs(fAv));
|
||||||
REQUIRE(fabs(fAv - tAt.InnerProduct(vt, ft)) < tol * fabs(fAv));
|
REQUIRE(fabs(fAv - tAt.InnerProduct(vt, ft)) < tol * fabs(fAv));
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
@@ -185,9 +185,9 @@ TEST_CASE("DoF Transformation Functions",
|
|||||||
|
|
||||||
double fAv = A.InnerProduct(v, f);
|
double fAv = A.InnerProduct(v, f);
|
||||||
|
|
||||||
REQUIRE(fabs(fAv - nAn.InnerProduct(v , f )) < tol * fabs(fAv));
|
REQUIRE(fabs(fAv - nAn.InnerProduct(v, f )) < tol * fabs(fAv));
|
||||||
REQUIRE(fabs(fAv - At.InnerProduct(vt, f )) < tol * fabs(fAv));
|
REQUIRE(fabs(fAv - At.InnerProduct(vt, f )) < tol * fabs(fAv));
|
||||||
REQUIRE(fabs(fAv - tA.InnerProduct(v , ft)) < tol * fabs(fAv));
|
REQUIRE(fabs(fAv - tA.InnerProduct(v, ft)) < tol * fabs(fAv));
|
||||||
REQUIRE(fabs(fAv - tAt.InnerProduct(vt, ft)) < tol * fabs(fAv));
|
REQUIRE(fabs(fAv - tAt.InnerProduct(vt, ft)) < tol * fabs(fAv));
|
||||||
}
|
}
|
||||||
SECTION("TransformDual")
|
SECTION("TransformDual")
|
||||||
@@ -217,9 +217,9 @@ TEST_CASE("DoF Transformation Functions",
|
|||||||
|
|
||||||
double uAv = A.InnerProduct(v, u);
|
double uAv = A.InnerProduct(v, u);
|
||||||
|
|
||||||
REQUIRE(fabs(uAv - nAn.InnerProduct(v , u )) < tol * fabs(uAv));
|
REQUIRE(fabs(uAv - nAn.InnerProduct(v, u )) < tol * fabs(uAv));
|
||||||
REQUIRE(fabs(uAv - At.InnerProduct(vt, u )) < tol * fabs(uAv));
|
REQUIRE(fabs(uAv - At.InnerProduct(vt, u )) < tol * fabs(uAv));
|
||||||
REQUIRE(fabs(uAv - tA.InnerProduct(v , ut)) < tol * fabs(uAv));
|
REQUIRE(fabs(uAv - tA.InnerProduct(v, ut)) < tol * fabs(uAv));
|
||||||
REQUIRE(fabs(uAv - tAt.InnerProduct(vt, ut)) < tol * fabs(uAv));
|
REQUIRE(fabs(uAv - tAt.InnerProduct(vt, ut)) < tol * fabs(uAv));
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|||||||
Reference in New Issue
Block a user