MFEM  v4.2.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
staticcond.cpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2020, 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 "staticcond.hpp"
13 
14 namespace mfem
15 {
16 
18  : fes(fespace)
19 {
20  tr_fec = fespace->FEColl()->GetTraceCollection();
21  int vdim = fes->GetVDim();
22  int ordering = fes->GetOrdering();
23 #ifndef MFEM_USE_MPI
24  tr_fes = new FiniteElementSpace(fes->GetMesh(), tr_fec, vdim, ordering);
25 #else
26  pfes = dynamic_cast<ParFiniteElementSpace*>(fes);
27  if (!pfes)
28  {
29  tr_fes = new FiniteElementSpace(fes->GetMesh(), tr_fec, vdim, ordering);
30  tr_pfes = NULL;
31  }
32  else
33  {
34  tr_pfes = new ParFiniteElementSpace(pfes->GetParMesh(), tr_fec, vdim,
35  ordering);
36  tr_fes = tr_pfes;
37  }
40 #endif
41  S = S_e = NULL;
42  symm = false;
43  A_data.Reset();
44  A_ipiv.Reset();
45 
46  Array<int> vdofs;
47  const int NE = fes->GetNE();
48  elem_pdof.MakeI(NE);
49  for (int i = 0; i < NE; i++)
50  {
51  const int npd = vdim*fes->GetNumElementInteriorDofs(i);
52  elem_pdof.AddColumnsInRow(i, npd);
53  }
54  elem_pdof.MakeJ();
55  for (int i = 0; i < NE; i++)
56  {
57  fes->GetElementVDofs(i, vdofs);
58  const int nsd = vdofs.Size()/vdim;
59  const int nspd = fes->GetNumElementInteriorDofs(i);
60  const int *dofs = vdofs.GetData();
61  for (int vd = 0; vd < vdim; vd++)
62  {
63 #ifdef MFEM_DEBUG
64  for (int j = 0; j < nspd; j++)
65  {
66  MFEM_ASSERT(dofs[nsd-nspd+j] >= 0, "");
67  }
68 #endif
69  elem_pdof.AddConnections(i, dofs+nsd-nspd, nspd);
70  dofs += nsd;
71  }
72  }
73  elem_pdof.ShiftUpI();
74  // Set the number of private dofs.
75  npdofs = elem_pdof.Size_of_connections();
76  MFEM_ASSERT(fes->GetVSize() == tr_fes->GetVSize() + npdofs,
77  "incompatible volume and trace FE spaces");
78  // Initialize the map rdof_edof.
79  rdof_edof.SetSize(tr_fes->GetVSize());
80  Array<int> rvdofs;
81  for (int i = 0; i < NE; i++)
82  {
83  fes->GetElementVDofs(i, vdofs);
84  tr_fes->GetElementVDofs(i, rvdofs);
85  const int nsd = vdofs.Size()/vdim;
86  const int nsrd = rvdofs.Size()/vdim;
87  for (int vd = 0; vd < vdim; vd++)
88  {
89  for (int j = 0; j < nsrd; j++)
90  {
91  int rvdof = rvdofs[j+nsrd*vd];
92  int vdof = vdofs[j+nsd*vd];
93  if (rvdof < 0)
94  {
95  rvdof = -1-rvdof;
96  vdof = -1-vdof;
97  }
98  MFEM_ASSERT(vdof >= 0, "incompatible volume and trace FE spaces");
99  rdof_edof[rvdof] = vdof;
100  }
101  }
102  }
103 }
104 
106 {
107 #ifdef MFEM_USE_MPI
108  // pS, pS_e are automatically destroyed
109 #endif
110  delete S_e;
111  delete S;
112  A_data.Delete();
113  A_ipiv.Delete();
114  delete tr_fes;
115  delete tr_fec;
116 }
117 
119 {
120  if (!Parallel())
121  {
122  return (tr_fes->GetTrueVSize() < fes->GetTrueVSize());
123  }
124  else
125  {
126 #ifdef MFEM_USE_MPI
127  return (tr_pfes->GlobalTrueVSize() < pfes->GlobalTrueVSize());
128 #else
129  return false; // avoid compiler warning
130 #endif
131  }
132 }
133 
134 void StaticCondensation::Init(bool symmetric, bool block_diagonal)
135 {
136  const int NE = fes->GetNE();
137  // symm = symmetric; // TODO: handle the symmetric case
138  A_offsets.SetSize(NE+1);
139  A_ipiv_offsets.SetSize(NE+1);
140  A_offsets[0] = A_ipiv_offsets[0] = 0;
141  Array<int> rvdofs;
142  for (int i = 0; i < NE; i++)
143  {
144  tr_fes->GetElementVDofs(i, rvdofs);
145  const int ned = rvdofs.Size();
146  const int npd = elem_pdof.RowSize(i);
147  A_offsets[i+1] = A_offsets[i] + npd*(npd + (symm ? 1 : 2)*ned);
148  A_ipiv_offsets[i+1] = A_ipiv_offsets[i] + npd;
149  }
150  A_data = Memory<double>(A_offsets[NE]);
151  A_ipiv = Memory<int>(A_ipiv_offsets[NE]);
152  const int nedofs = tr_fes->GetVSize();
153  if (fes->GetVDim() == 1)
154  {
155  // The sparsity pattern of S is given by the map rdof->elem->rdof
156  Table rdof_rdof;
157  {
158  Table elem_rdof, rdof_elem;
159  elem_rdof.MakeI(NE);
160  for (int i = 0; i < NE; i++)
161  {
162  tr_fes->GetElementVDofs(i, rvdofs);
163  elem_rdof.AddColumnsInRow(i, rvdofs.Size());
164  }
165  elem_rdof.MakeJ();
166  for (int i = 0; i < NE; i++)
167  {
168  tr_fes->GetElementVDofs(i, rvdofs);
170  elem_rdof.AddConnections(i, rvdofs.GetData(), rvdofs.Size());
171  }
172  elem_rdof.ShiftUpI();
173  Transpose(elem_rdof, rdof_elem, nedofs);
174  mfem::Mult(rdof_elem, elem_rdof, rdof_rdof);
175  }
176  S = new SparseMatrix(rdof_rdof.GetI(), rdof_rdof.GetJ(), NULL,
177  nedofs, nedofs, true, false, false);
178  rdof_rdof.LoseData();
179  }
180  else
181  {
182  // For a block diagonal vector bilinear form, the sparsity of
183  // rdof->elem->rdof is overkill, so we use dynamically allocated
184  // sparsity pattern.
185  S = new SparseMatrix(nedofs);
186  }
187 }
188 
190 {
191  Array<int> rvdofs;
192  tr_fes->GetElementVDofs(el, rvdofs);
193  const int vdim = fes->GetVDim();
194  const int nvpd = elem_pdof.RowSize(el);
195  const int nved = rvdofs.Size();
196  DenseMatrix A_pp(A_data + A_offsets[el], nvpd, nvpd);
197  DenseMatrix A_pe(A_pp.Data() + nvpd*nvpd, nvpd, nved);
198  DenseMatrix A_ep;
199  if (symm) { A_ep.SetSize(nved, nvpd); }
200  else { A_ep.UseExternalData(A_pe.Data() + nvpd*nved, nved, nvpd); }
201  DenseMatrix A_ee(nved, nved);
202 
203  const int npd = nvpd/vdim;
204  const int ned = nved/vdim;
205  const int nd = npd + ned;
206  // Copy the blocks from elmat to A_xx
207  for (int i = 0; i < vdim; i++)
208  {
209  for (int j = 0; j < vdim; j++)
210  {
211  A_pp.CopyMN(elmat, npd, npd, i*nd+ned, j*nd+ned, i*npd, j*npd);
212  A_pe.CopyMN(elmat, npd, ned, i*nd+ned, j*nd, i*npd, j*ned);
213  A_ep.CopyMN(elmat, ned, npd, i*nd, j*nd+ned, i*ned, j*npd);
214  A_ee.CopyMN(elmat, ned, ned, i*nd, j*nd, i*ned, j*ned);
215  }
216  }
217  // Compute the Schur complement
218  LUFactors lu(A_pp.Data(), A_ipiv + A_ipiv_offsets[el]);
219  lu.Factor(nvpd);
220  lu.BlockFactor(nvpd, nved, A_pe.Data(), A_ep.Data(), A_ee.Data());
221 
222  // Assemble the Schur complement
223  const int skip_zeros = 0;
224  S->AddSubMatrix(rvdofs, rvdofs, A_ee, skip_zeros);
225 }
226 
228 {
229  Array<int> rvdofs;
230  tr_fes->GetBdrElementVDofs(el, rvdofs);
231  const int skip_zeros = 0;
232  S->AddSubMatrix(rvdofs, rvdofs, elmat, skip_zeros);
233 }
234 
236 {
237  const int skip_zeros = 0;
238  if (!Parallel())
239  {
240  S->Finalize(skip_zeros);
241  if (S_e) { S_e->Finalize(skip_zeros); }
242  const SparseMatrix *cP = tr_fes->GetConformingProlongation();
243  if (cP)
244  {
245  if (S->Height() != cP->Width())
246  {
247  SparseMatrix *cS = mfem::RAP(*cP, *S, *cP);
248  delete S;
249  S = cS;
250  }
251  if (S_e && S_e->Height() != cP->Width())
252  {
253  SparseMatrix *cS_e = mfem::RAP(*cP, *S_e, *cP);
254  delete S_e;
255  S = cS_e;
256  }
257  }
258  }
259  else // parallel
260  {
261 #ifdef MFEM_USE_MPI
262  if (!S) { return; } // already finalized
263  S->Finalize(skip_zeros);
264  if (S_e) { S_e->Finalize(skip_zeros); }
265  OperatorHandle dS(pS.Type()), pP(pS.Type());
266  dS.MakeSquareBlockDiag(tr_pfes->GetComm(), tr_pfes->GlobalVSize(),
267  tr_pfes->GetDofOffsets(), S);
268  // TODO - construct Dof_TrueDof_Matrix directly in the pS format
269  pP.ConvertFrom(tr_pfes->Dof_TrueDof_Matrix());
270  pS.MakePtAP(dS, pP);
271  dS.Clear();
272  delete S;
273  S = NULL;
274  if (S_e)
275  {
276  OperatorHandle dS_e(pS_e.Type());
277  dS_e.MakeSquareBlockDiag(tr_pfes->GetComm(), tr_pfes->GlobalVSize(),
278  tr_pfes->GetDofOffsets(), S_e);
279  pS_e.MakePtAP(dS_e, pP);
280  dS_e.Clear();
281  delete S_e;
282  S_e = NULL;
283  }
284 #endif
285  }
286 }
287 
289  const Array<int> &ess_rtdof_list, Matrix::DiagonalPolicy dpolicy)
290 {
291  if (!Parallel() || S) // not parallel or not finalized
292  {
293  if (S_e == NULL)
294  {
295  S_e = new SparseMatrix(S->Height());
296  }
297  for (int i = 0; i < ess_rtdof_list.Size(); i++)
298  {
299  S->EliminateRowCol(ess_rtdof_list[i], *S_e, dpolicy);
300  }
301  }
302  else // parallel and finalized
303  {
304 #ifdef MFEM_USE_MPI
305  MFEM_ASSERT(pS_e.Ptr() == NULL, "essential b.c. already eliminated");
306  pS_e.EliminateRowsCols(pS, ess_rtdof_list);
307 #endif
308  }
309 }
310 
311 void StaticCondensation::ReduceRHS(const Vector &b, Vector &sc_b) const
312 {
313  // sc_b = b_e - A_ep A_pp_inv b_p
314 
315  MFEM_ASSERT(b.Size() == fes->GetVSize(), "'b' has incorrect size");
316 
317  const int NE = fes->GetNE();
318  const int nedofs = tr_fes->GetVSize();
319  const SparseMatrix *tr_cP = NULL;
320  Vector b_r;
321  if (!Parallel() && !(tr_cP = tr_fes->GetConformingProlongation()))
322  {
323  sc_b.SetSize(nedofs);
324  b_r.SetDataAndSize(sc_b.GetData(), sc_b.Size());
325  }
326  else
327  {
328  b_r.SetSize(nedofs);
329  }
330  for (int i = 0; i < nedofs; i++)
331  {
332  b_r(i) = b(rdof_edof[i]);
333  }
334 
335  DenseMatrix U_pe, L_ep;
336  Vector b_p, b_ep;
337  Array<int> rvdofs;
338  for (int i = 0; i < NE; i++)
339  {
340  tr_fes->GetElementVDofs(i, rvdofs);
341  const int ned = rvdofs.Size();
342  const int *rd = rvdofs.GetData();
343  const int npd = elem_pdof.RowSize(i);
344  const int *pd = elem_pdof.GetRow(i);
345  b_p.SetSize(npd);
346  b_ep.SetSize(ned);
347  for (int j = 0; j < npd; j++)
348  {
349  b_p(j) = b(pd[j]);
350  }
351 
352  LUFactors lu(const_cast<double*>((const double*)A_data) + A_offsets[i],
353  const_cast<int*>((const int*)A_ipiv) + A_ipiv_offsets[i]);
354  lu.LSolve(npd, 1, b_p);
355 
356  if (symm)
357  {
358  // TODO: handle the symmetric case correctly.
359  U_pe.UseExternalData(lu.data + npd*npd, npd, ned);
360  U_pe.MultTranspose(b_p, b_ep);
361  }
362  else
363  {
364  L_ep.UseExternalData(lu.data + npd*(npd+ned), ned, npd);
365  L_ep.Mult(b_p, b_ep);
366  }
367  for (int j = 0; j < ned; j++)
368  {
369  if (rd[j] >= 0) { b_r(rd[j]) -= b_ep(j); }
370  else { b_r(-1-rd[j]) += b_ep(j); }
371  }
372  }
373  if (!Parallel())
374  {
375  if (tr_cP)
376  {
377  sc_b.SetSize(tr_cP->Width());
378  tr_cP->MultTranspose(b_r, sc_b);
379  }
380  }
381  else
382  {
383 #ifdef MFEM_USE_MPI
384  const Operator *tr_P = tr_pfes->GetProlongationMatrix();
385  sc_b.SetSize(tr_P->Width());
386  tr_P->MultTranspose(b_r, sc_b);
387 #endif
388  }
389 }
390 
391 void StaticCondensation::ReduceSolution(const Vector &sol, Vector &sc_sol) const
392 {
393  MFEM_ASSERT(sol.Size() == fes->GetVSize(), "'sol' has incorrect size");
394 
395  const int nedofs = tr_fes->GetVSize();
396  const SparseMatrix *tr_R = tr_fes->GetRestrictionMatrix();
397  Vector sol_r;
398  if (!tr_R)
399  {
400  sc_sol.SetSize(nedofs);
401  sol_r.SetDataAndSize(sc_sol.GetData(), sc_sol.Size());
402  }
403  else
404  {
405  sol_r.SetSize(nedofs);
406  }
407  for (int i = 0; i < nedofs; i++)
408  {
409  sol_r(i) = sol(rdof_edof[i]);
410  }
411  if (tr_R)
412  {
413  sc_sol.SetSize(tr_R->Height());
414  tr_R->Mult(sol_r, sc_sol);
415  }
416 }
417 
419  Vector &B, int copy_interior) const
420 {
421  ReduceRHS(b, B);
422  ReduceSolution(x, X);
423  if (!Parallel())
424  {
425  S_e->AddMult(X, B, -1.);
426  S->PartMult(ess_rtdof_list, X, B);
427  }
428  else
429  {
430 #ifdef MFEM_USE_MPI
431  MFEM_ASSERT(pS.Type() == pS_e.Type(), "type id mismatch");
432  pS.EliminateBC(pS_e, ess_rtdof_list, X, B);
433 #endif
434  }
435  if (!copy_interior)
436  {
437  X.SetSubVectorComplement(ess_rtdof_list, 0.0);
438  }
439 }
440 
442  const Array<int> &ess_tdof_marker, Array<int> &ess_rtdof_marker) const
443 {
444  const int nedofs = tr_fes->GetVSize();
445  const SparseMatrix *R = fes->GetRestrictionMatrix();
446  Array<int> ess_dof_marker;
447  if (!R)
448  {
449  ess_dof_marker.MakeRef(ess_tdof_marker);
450  }
451  else
452  {
453  ess_dof_marker.SetSize(fes->GetVSize());
454  R->BooleanMultTranspose(ess_tdof_marker, ess_dof_marker);
455  }
456  const SparseMatrix *tr_R = tr_fes->GetRestrictionMatrix();
457  Array<int> ess_rdof_marker;
458  if (!tr_R)
459  {
460  ess_rtdof_marker.SetSize(nedofs);
461  ess_rdof_marker.MakeRef(ess_rtdof_marker);
462  }
463  else
464  {
465  ess_rdof_marker.SetSize(nedofs);
466  }
467  for (int i = 0; i < nedofs; i++)
468  {
469  ess_rdof_marker[i] = ess_dof_marker[rdof_edof[i]];
470  }
471  if (tr_R)
472  {
473  ess_rtdof_marker.SetSize(tr_R->Height());
474  tr_R->BooleanMult(ess_rdof_marker, ess_rtdof_marker);
475  }
476 }
477 
479  const Vector &b, const Vector &sc_sol, Vector &sol) const
480 {
481  // sol_e = sc_sol
482  // sol_p = A_pp_inv (b_p - A_pe sc_sol)
483 
484  MFEM_ASSERT(b.Size() == fes->GetVSize(), "'b' has incorrect size");
485 
486  const int nedofs = tr_fes->GetVSize();
487  Vector sol_r;
488  if (!Parallel())
489  {
490  const SparseMatrix *tr_cP = tr_fes->GetConformingProlongation();
491  if (!tr_cP)
492  {
493  sol_r.SetDataAndSize(sc_sol.GetData(), sc_sol.Size());
494  }
495  else
496  {
497  sol_r.SetSize(nedofs);
498  tr_cP->Mult(sc_sol, sol_r);
499  }
500  }
501  else
502  {
503 #ifdef MFEM_USE_MPI
504  sol_r.SetSize(nedofs);
505  tr_pfes->GetProlongationMatrix()->Mult(sc_sol, sol_r);
506 #endif
507  }
508  sol.SetSize(nedofs+npdofs);
509  for (int i = 0; i < nedofs; i++)
510  {
511  sol(rdof_edof[i]) = sol_r(i);
512  }
513  const int NE = fes->GetNE();
514  Vector b_p, s_e;
515  Array<int> rvdofs;
516  for (int i = 0; i < NE; i++)
517  {
518  tr_fes->GetElementVDofs(i, rvdofs);
519  const int ned = rvdofs.Size();
520  const int npd = elem_pdof.RowSize(i);
521  const int *pd = elem_pdof.GetRow(i);
522  b_p.SetSize(npd);
523 
524  for (int j = 0; j < npd; j++)
525  {
526  b_p(j) = b(pd[j]);
527  }
528  sol_r.GetSubVector(rvdofs, s_e);
529 
530  LUFactors lu(const_cast<double*>((const double*)A_data) + A_offsets[i],
531  const_cast<int*>((const int*)A_ipiv) + A_ipiv_offsets[i]);
532  lu.LSolve(npd, 1, b_p);
533  lu.BlockBackSolve(npd, ned, 1, lu.data + npd*npd, s_e, b_p);
534 
535  for (int j = 0; j < npd; j++)
536  {
537  sol(pd[j]) = b_p(j);
538  }
539  }
540 }
541 
542 }
Ordering::Type GetOrdering() const
Return the ordering method.
Definition: fespace.hpp:412
int Size() const
Return the logical size of the array.
Definition: array.hpp:124
int GetVSize() const
Return the number of vector dofs, i.e. GetNDofs() x GetVDim().
Definition: fespace.hpp:400
int * GetJ()
Definition: table.hpp:114
virtual FiniteElementCollection * GetTraceCollection() const
Definition: fe_coll.cpp:41
const SparseMatrix * GetConformingProlongation() const
The returned SparseMatrix is owned by the FiniteElementSpace.
Definition: fespace.cpp:874
~StaticCondensation()
Destroy a StaticCondensation object.
Definition: staticcond.cpp:105
void AddColumnsInRow(int r, int ncol)
Definition: table.hpp:78
bool ReducesTrueVSize() const
Definition: staticcond.cpp:118
HypreParMatrix * RAP(const HypreParMatrix *A, const HypreParMatrix *P)
Returns the matrix P^t * A * P.
Definition: hypre.cpp:1640
void ReduceRHS(const Vector &b, Vector &sc_b) const
Definition: staticcond.cpp:311
void MakeI(int nrows)
Next 7 methods are used together with the default constructor.
Definition: table.cpp:85
HYPRE_Int * GetDofOffsets() const
Definition: pfespace.hpp:247
virtual void Finalize(int skip_zeros=1)
Finalize the matrix initialization, switching the storage format from LIL to CSR. ...
Definition: sparsemat.hpp:460
virtual const Operator * GetProlongationMatrix() const
The returned Operator is owned by the FiniteElementSpace.
Definition: pfespace.cpp:895
void SetSize(int s)
Resize the vector to size s.
Definition: vector.hpp:459
void GetElementVDofs(int i, Array< int > &vdofs) const
Returns indexes of degrees of freedom in array dofs for i&#39;th element.
Definition: fespace.cpp:173
void Delete()
Delete the owned pointers. The Memory is not reset by this method, i.e. it will, generally, not be Empty() after this call.
void BooleanMult(const Array< int > &x, Array< int > &y) const
y = A * x, treating all entries as booleans (zero=false, nonzero=true).
Definition: sparsemat.cpp:831
void Mult(const Table &A, const Table &B, Table &C)
C = A * B (as boolean matrices)
Definition: table.cpp:476
int Width() const
Get the width (size of input) of the Operator. Synonym with NumCols().
Definition: operator.hpp:71
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
Operator::Type Type() const
Get the currently set operator type id.
Definition: handle.hpp:91
void BlockBackSolve(int m, int n, int r, const double *U12, const double *X2, double *Y1) const
Definition: densemat.cpp:3208
void EliminateBC(const OperatorHandle &A_e, const Array< int > &ess_dof_list, const Vector &X, Vector &B) const
Eliminate essential dofs from the solution X into the r.h.s. B.
Definition: handle.cpp:340
void GetSubVector(const Array< int > &dofs, Vector &elemvect) const
Extract entries listed in dofs to the output Vector elemvect.
Definition: vector.cpp:495
T * GetData()
Returns the data.
Definition: array.hpp:98
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
void AddMult(const Vector &x, Vector &y, const double a=1.0) const
y += A * x (default) or y += a * A * x
Definition: sparsemat.cpp:601
int Size() const
Returns the size of the vector.
Definition: vector.hpp:160
void GetRow(int i, Array< int > &row) const
Return row i in array row (the Table must be finalized)
Definition: table.cpp:191
int GetNumElementInteriorDofs(int i) const
Definition: fespace.hpp:496
void ReduceSystem(Vector &x, Vector &b, Vector &X, Vector &B, int copy_interior=0) const
Set the reduced solution X and r.h.s B vectors from the full linear system solution x and r...
Definition: staticcond.cpp:418
void EliminateReducedTrueDofs(const Array< int > &ess_rtdof_list, Matrix::DiagonalPolicy dpolicy)
Eliminate the given reduced true dofs from the Schur complement matrix S.
Definition: staticcond.cpp:288
void LoseData()
Call this if data has been stolen.
Definition: table.hpp:154
virtual void Mult(const Vector &x, Vector &y) const =0
Operator application: y=A(x).
Abstract parallel finite element space.
Definition: pfespace.hpp:28
StaticCondensation(FiniteElementSpace *fespace)
Construct a StaticCondensation object.
Definition: staticcond.cpp:17
void ComputeSolution(const Vector &b, const Vector &sc_sol, Vector &sol) const
Definition: staticcond.cpp:478
double * GetData() const
Return a pointer to the beginning of the Vector data.
Definition: vector.hpp:169
int Size_of_connections() const
Definition: table.hpp:98
void BooleanMultTranspose(const Array< int > &x, Array< int > &y) const
y = At * x, treating all entries as booleans (zero=false, nonzero=true).
Definition: sparsemat.cpp:861
void AddConnections(int r, const int *c, int nc)
Definition: table.cpp:108
virtual void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
Definition: operator.hpp:92
int GetNE() const
Returns number of elements in the mesh.
Definition: fespace.hpp:427
void MultTranspose(const double *x, double *y) const
Multiply a vector with the transpose matrix.
Definition: densemat.cpp:203
MPI_Comm GetComm() const
Definition: pfespace.hpp:239
Data type sparse matrix.
Definition: sparsemat.hpp:46
int Height() const
Get the height (size of output) of the Operator. Synonym with NumRows().
Definition: operator.hpp:65
Mesh * GetMesh() const
Returns the mesh.
Definition: fespace.hpp:314
double b
Definition: lissajous.cpp:42
void EliminateRowsCols(OperatorHandle &A, const Array< int > &ess_dof_list)
Reset the OperatorHandle to be the eliminated part of A after elimination of the essential dofs ess_d...
Definition: handle.cpp:252
HypreParMatrix * Dof_TrueDof_Matrix() const
The true dof-to-dof interpolation matrix.
Definition: pfespace.hpp:281
HYPRE_Int GlobalTrueVSize() const
Definition: pfespace.hpp:251
void LSolve(int m, int n, double *X) const
Definition: densemat.cpp:2977
void Transpose(const Table &A, Table &At, int _ncols_A)
Transpose a Table.
Definition: table.cpp:414
virtual int GetTrueVSize() const
Return the number of vector true (conforming) dofs.
Definition: fespace.hpp:403
int GetVDim() const
Returns vector dimension.
Definition: fespace.hpp:389
void Reset()
Reset the memory to be empty, ensuring that Delete() will be a no-op.
void AddSubMatrix(const Array< int > &rows, const Array< int > &cols, const DenseMatrix &subm, int skip_zeros=1)
Definition: sparsemat.cpp:2436
Operator * Ptr() const
Access the underlying Operator pointer.
Definition: handle.hpp:82
void SetSubVectorComplement(const Array< int > &dofs, const double val)
Set all vector entries NOT in the dofs Array to the given val.
Definition: vector.cpp:656
HYPRE_Int GlobalVSize() const
Definition: pfespace.hpp:249
double * Data() const
Returns the matrix data array.
Definition: densemat.hpp:96
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:87
virtual void Mult(const Vector &x, Vector &y) const
Matrix vector multiplication.
Definition: sparsemat.cpp:594
void Finalize()
Finalize the construction of the Schur complement matrix.
Definition: staticcond.cpp:235
void SetSize(int nsize)
Change the logical size of the array, keep existing entries.
Definition: array.hpp:654
void EliminateRowCol(int rc, const double sol, Vector &rhs, DiagonalPolicy dpolicy=DIAG_ONE)
Eliminate row rc and column rc and modify the rhs using sol.
Definition: sparsemat.cpp:1624
void SetDataAndSize(double *d, int s)
Set the Vector data and size.
Definition: vector.hpp:128
void AssembleBdrMatrix(int el, const DenseMatrix &elmat)
Definition: staticcond.cpp:227
void ShiftUpI()
Definition: table.cpp:119
void MakeSquareBlockDiag(MPI_Comm comm, HYPRE_Int glob_size, HYPRE_Int *row_starts, SparseMatrix *diag)
Reset the OperatorHandle to hold a parallel square block-diagonal matrix using the currently set type...
Definition: handle.cpp:60
void AssembleMatrix(int el, const DenseMatrix &elmat)
Definition: staticcond.cpp:189
DiagonalPolicy
Defines operator diagonal policy upon elimination of rows and/or columns.
Definition: operator.hpp:45
void MakeJ()
Definition: table.cpp:95
bool Factor(int m, double TOL=0.0)
Compute the LU factorization of the current matrix.
Definition: densemat.cpp:2865
void MultTranspose(const Vector &x, Vector &y) const
Multiply a vector with the transposed matrix. y = At * x.
Definition: sparsemat.cpp:721
void MakeRef(T *, int)
Make this Array a reference to a pointer.
Definition: array.hpp:839
Vector data type.
Definition: vector.hpp:51
const FiniteElementCollection * FEColl() const
Definition: fespace.hpp:414
void Mult(const double *x, double *y) const
Matrix vector multiplication.
Definition: densemat.cpp:175
ID for class HypreParMatrix.
Definition: operator.hpp:241
int * GetI()
Definition: table.hpp:113
void CopyMN(const DenseMatrix &A, int m, int n, int Aro, int Aco)
Copy the m x n submatrix of A at row/col offsets Aro/Aco to *this.
Definition: densemat.cpp:1536
void ReduceSolution(const Vector &sol, Vector &sc_sol) const
Definition: staticcond.cpp:391
void UseExternalData(double *d, int h, int w)
Change the data array and the size of the DenseMatrix.
Definition: densemat.hpp:65
int RowSize(int i) const
Definition: table.hpp:108
void GetBdrElementVDofs(int i, Array< int > &vdofs) const
Returns indexes of degrees of freedom for i&#39;th boundary element.
Definition: fespace.cpp:179
void SetSize(int s)
Change the size of the DenseMatrix to s x s.
Definition: densemat.hpp:90
Abstract operator.
Definition: operator.hpp:24
virtual const SparseMatrix * GetRestrictionMatrix() const
The returned SparseMatrix is owned by the FiniteElementSpace.
Definition: fespace.hpp:334
void Init(bool symmetric, bool block_diagonal)
Definition: staticcond.cpp:134
void ConvertMarkerToReducedTrueDofs(const Array< int > &ess_tdof_marker, Array< int > &ess_rtdof_marker) const
Definition: staticcond.cpp:441
void MakePtAP(OperatorHandle &A, OperatorHandle &P)
Reset the OperatorHandle to hold the product P^t A P.
Definition: handle.cpp:123
static void AdjustVDofs(Array< int > &vdofs)
Definition: fespace.cpp:160
void SetType(Operator::Type tid)
Invoke Clear() and set a new type id.
Definition: handle.hpp:124
void PartMult(const Array< int > &rows, const Vector &x, Vector &y) const
Definition: sparsemat.cpp:787
double * data
Definition: densemat.hpp:515