MFEM  v4.5.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
bilininteg_transpose_ea.cpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2022, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #include "../general/forall.hpp"
13 #include "bilininteg.hpp"
14 
15 namespace mfem
16 {
17 
19  Vector &ea_data, const bool add)
20 {
21  if (add)
22  {
23  Vector ea_data_tmp(ea_data.Size());
24  bfi->AssembleEA(fes, ea_data_tmp, false);
25  const int ne = fes.GetNE();
26  if (ne == 0) { return; }
27  const int dofs = fes.GetFE(0)->GetDof();
28  auto A = Reshape(ea_data_tmp.Read(), dofs, dofs, ne);
29  auto AT = Reshape(ea_data.ReadWrite(), dofs, dofs, ne);
30  MFEM_FORALL(e, ne,
31  {
32  for (int i = 0; i < dofs; i++)
33  {
34  for (int j = 0; j < dofs; j++)
35  {
36  const double a = A(i, j, e);
37  AT(j, i, e) += a;
38  }
39  }
40  });
41  }
42  else
43  {
44  bfi->AssembleEA(fes, ea_data, false);
45  const int ne = fes.GetNE();
46  if (ne == 0) { return; }
47  const int dofs = fes.GetFE(0)->GetDof();
48  auto A = Reshape(ea_data.ReadWrite(), dofs, dofs, ne);
49  MFEM_FORALL(e, ne,
50  {
51  for (int i = 0; i < dofs; i++)
52  {
53  for (int j = i+1; j < dofs; j++)
54  {
55  const double aij = A(i, j, e);
56  const double aji = A(j, i, e);
57  A(j, i, e) = aij;
58  A(i, j, e) = aji;
59  }
60  }
61  });
62  }
63 }
64 
66  Vector &ea_data_int,
67  Vector &ea_data_ext,
68  const bool add)
69 {
70  const int nf = fes.GetNFbyType(FaceType::Interior);
71  if (nf == 0) { return; }
72  if (add)
73  {
74  Vector ea_data_int_tmp(ea_data_int.Size());
75  Vector ea_data_ext_tmp(ea_data_ext.Size());
76  bfi->AssembleEAInteriorFaces(fes, ea_data_int_tmp, ea_data_ext_tmp, false);
77  const int faceDofs = fes.GetTraceElement(0,
78  fes.GetMesh()->GetFaceBaseGeometry(0))->GetDof();
79  auto A_int = Reshape(ea_data_int_tmp.Read(), faceDofs, faceDofs, 2, nf);
80  auto A_ext = Reshape(ea_data_ext_tmp.Read(), faceDofs, faceDofs, 2, nf);
81  auto AT_int = Reshape(ea_data_int.ReadWrite(), faceDofs, faceDofs, 2, nf);
82  auto AT_ext = Reshape(ea_data_ext.ReadWrite(), faceDofs, faceDofs, 2, nf);
83  MFEM_FORALL(f, nf,
84  {
85  for (int i = 0; i < faceDofs; i++)
86  {
87  for (int j = 0; j < faceDofs; j++)
88  {
89  const double a_int0 = A_int(i, j, 0, f);
90  const double a_int1 = A_int(i, j, 1, f);
91  const double a_ext0 = A_ext(i, j, 0, f);
92  const double a_ext1 = A_ext(i, j, 1, f);
93  AT_int(j, i, 0, f) += a_int0;
94  AT_int(j, i, 1, f) += a_int1;
95  AT_ext(j, i, 0, f) += a_ext1;
96  AT_ext(j, i, 1, f) += a_ext0;
97  }
98  }
99  });
100  }
101  else
102  {
103  bfi->AssembleEAInteriorFaces(fes, ea_data_int, ea_data_ext, false);
104  const int faceDofs = fes.GetTraceElement(0,
105  fes.GetMesh()->GetFaceBaseGeometry(0))->GetDof();
106  auto A_int = Reshape(ea_data_int.ReadWrite(), faceDofs, faceDofs, 2, nf);
107  auto A_ext = Reshape(ea_data_ext.ReadWrite(), faceDofs, faceDofs, 2, nf);
108  MFEM_FORALL(f, nf,
109  {
110  for (int i = 0; i < faceDofs; i++)
111  {
112  for (int j = i+1; j < faceDofs; j++)
113  {
114  const double aij_int0 = A_int(i, j, 0, f);
115  const double aij_int1 = A_int(i, j, 1, f);
116  const double aji_int0 = A_int(j, i, 0, f);
117  const double aji_int1 = A_int(j, i, 1, f);
118  A_int(j, i, 0, f) = aij_int0;
119  A_int(j, i, 1, f) = aij_int1;
120  A_int(i, j, 0, f) = aji_int0;
121  A_int(i, j, 1, f) = aji_int1;
122  }
123  }
124  for (int i = 0; i < faceDofs; i++)
125  {
126  for (int j = 0; j < faceDofs; j++)
127  {
128  const double aij_ext0 = A_ext(i, j, 0, f);
129  const double aji_ext1 = A_ext(j, i, 1, f);
130  A_ext(j, i, 1, f) = aij_ext0;
131  A_ext(i, j, 0, f) = aji_ext1;
132  }
133  }
134  });
135  }
136 }
137 
139  Vector &ea_data_bdr,
140  const bool add)
141 {
142  const int nf = fes.GetNFbyType(FaceType::Boundary);
143  if (nf == 0) { return; }
144  if (add)
145  {
146  Vector ea_data_bdr_tmp(ea_data_bdr.Size());
147  bfi->AssembleEABoundaryFaces(fes, ea_data_bdr_tmp, false);
148  const int faceDofs = fes.GetTraceElement(0,
149  fes.GetMesh()->GetFaceBaseGeometry(0))->GetDof();
150  auto A_bdr = Reshape(ea_data_bdr_tmp.Read(), faceDofs, faceDofs, nf);
151  auto AT_bdr = Reshape(ea_data_bdr.ReadWrite(), faceDofs, faceDofs, nf);
152  MFEM_FORALL(f, nf,
153  {
154  for (int i = 0; i < faceDofs; i++)
155  {
156  for (int j = 0; j < faceDofs; j++)
157  {
158  const double a_bdr = A_bdr(i, j, f);
159  AT_bdr(j, i, f) += a_bdr;
160  }
161  }
162  });
163  }
164  else
165  {
166  bfi->AssembleEABoundaryFaces(fes, ea_data_bdr, false);
167  const int faceDofs = fes.GetTraceElement(0,
168  fes.GetMesh()->GetFaceBaseGeometry(0))->GetDof();
169  auto A_bdr = Reshape(ea_data_bdr.ReadWrite(), faceDofs, faceDofs, nf);
170  MFEM_FORALL(f, nf,
171  {
172  for (int i = 0; i < faceDofs; i++)
173  {
174  for (int j = i+1; j < faceDofs; j++)
175  {
176  const double aij_bdr = A_bdr(i, j, f);
177  const double aji_bdr = A_bdr(j, i, f);
178  A_bdr(j, i, f) = aij_bdr;
179  A_bdr(i, j, f) = aji_bdr;
180  }
181  }
182  });
183  }
184 }
185 
186 }
int GetNFbyType(FaceType type) const
Returns the number of faces according to the requested type.
Definition: fespace.hpp:631
virtual void AssembleEABoundaryFaces(const FiniteElementSpace &fes, Vector &ea_data_bdr, const bool add=true)
Definition: bilininteg.cpp:72
const FiniteElement * GetTraceElement(int i, Geometry::Type geom_type) const
Return the trace element from element &#39;i&#39; to the given &#39;geom_type&#39;.
Definition: fespace.cpp:3173
virtual void AssembleEABoundaryFaces(const FiniteElementSpace &fes, Vector &ea_data_bdr, const bool add)
virtual void AssembleEA(const FiniteElementSpace &fes, Vector &emat, const bool add)
Method defining element assembly.
int Size() const
Returns the size of the vector.
Definition: vector.hpp:200
Geometry::Type GetFaceBaseGeometry(int i) const
Definition: mesh.hpp:1066
void add(const Vector &v1, const Vector &v2, Vector &v)
Definition: vector.cpp:313
virtual void AssembleEA(const FiniteElementSpace &fes, Vector &emat, const bool add=true)
Method defining element assembly.
Definition: bilininteg.cpp:54
int GetNE() const
Returns number of elements in the mesh.
Definition: fespace.hpp:614
double f(const Vector &xvec)
Definition: lor_mms.hpp:32
Mesh * GetMesh() const
Returns the mesh.
Definition: fespace.hpp:441
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:96
int GetDof() const
Returns the number of degrees of freedom in the finite element.
Definition: fe_base.hpp:323
virtual void AssembleEAInteriorFaces(const FiniteElementSpace &fes, Vector &ea_data_int, Vector &ea_data_ext, const bool add)
double a
Definition: lissajous.cpp:41
virtual double * ReadWrite(bool on_dev=true)
Shortcut for mfem::ReadWrite(vec.GetMemory(), vec.Size(), on_dev).
Definition: vector.hpp:465
virtual const FiniteElement * GetFE(int i) const
Returns pointer to the FiniteElement in the FiniteElementCollection associated with i&#39;th element in t...
Definition: fespace.cpp:2781
virtual void AssembleEAInteriorFaces(const FiniteElementSpace &fes, Vector &ea_data_int, Vector &ea_data_ext, const bool add=true)
Definition: bilininteg.cpp:62
Vector data type.
Definition: vector.hpp:60
MFEM_HOST_DEVICE DeviceTensor< sizeof...(Dims), T > Reshape(T *ptr, Dims...dims)
Wrap a pointer as a DeviceTensor with automatically deduced template parameters.
Definition: dtensor.hpp:131