Re: [deal.II] Direct solver matrix.free

2021-10-21 Thread Daniel Arndt
Joss,

The MatrixFree framework relies on the operator to be only defined by its
action. Accessing individual elements (as direct solvers typically require)
falls out of this concept. That is to say that the MatrixFree framework
only works well with an iterative solver.

Best,
Daniel

Am Do., 21. Okt. 2021 um 13:45 Uhr schrieb Joss G. :

> Hi all,
>
> I was using the solver PETScWrappers::SparseDirectMUMPS for a sparse
> matrix approach. Now I am changing the code to a matrix-free implementation
> and I am not using PETSc anymore, instead, I use
> LinearAlegebra::distributed::Vector. I am not able to find an
> equivalent direct solver that I can use. Could anyone please point me what
> solver can i use?
>
> Thank you very much.
>
> --
> The deal.II project is located at http://www.dealii.org/
> For mailing list/forum options, see
> https://groups.google.com/d/forum/dealii?hl=en
> ---
> You received this message because you are subscribed to the Google Groups
> "deal.II User Group" group.
> To unsubscribe from this group and stop receiving emails from it, send an
> email to dealii+unsubscr...@googlegroups.com.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/dealii/90916560-faf4-4d5a-b93c-d684deeb1b6fn%40googlegroups.com
> 
> .
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/CAOYDWbKw3LMWN1dRxWP2m1tHYmfJvy8JJMoyRURUTtwyk7pH2A%40mail.gmail.com.


[deal.II] Direct solver matrix.free

2021-10-21 Thread Joss G.
Hi all, 

I was using the solver PETScWrappers::SparseDirectMUMPS for a sparse matrix 
approach. Now I am changing the code to a matrix-free implementation and I 
am not using PETSc anymore, instead, I use 
LinearAlegebra::distributed::Vector. I am not able to find an 
equivalent direct solver that I can use. Could anyone please point me what 
solver can i use?

Thank you very much.

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/90916560-faf4-4d5a-b93c-d684deeb1b6fn%40googlegroups.com.


[deal.II] deal.II Newsletter #187

2021-10-21 Thread 'Rene Gassmoeller' via deal.II User Group
Hello everyone!

This is deal.II newsletter #187.
It automatically reports recently merged features and discussions about the 
deal.II finite element library.


## Below you find a list of recently proposed or merged features:

#12856: Fix typo (proposed by gfcas) https://github.com/dealii/dealii/pull/12856

#12855: fix step-49 geo file (proposed by tjhei) 
https://github.com/dealii/dealii/pull/12855

#12853: Implement convert_generalized_support_point_values_to_dof_values for 
wedges/pyramids (proposed by peterrum; merged) 
https://github.com/dealii/dealii/pull/12853

#12852: No more magic number Pi. (proposed by marcfehling; merged) 
https://github.com/dealii/dealii/pull/12852

#12851: DataOut: add wedge/pyramid test (proposed by peterrum) 
https://github.com/dealii/dealii/pull/12851

#12850: PreconditionChebyshev: Update documentation (proposed by kronbichler; 
merged) https://github.com/dealii/dealii/pull/12850

#12849: Reduce duplication of information. (proposed by bangerth; merged) 
https://github.com/dealii/dealii/pull/12849

#12848: Systematize construction of vertices for wedges. (proposed by bangerth; 
merged) https://github.com/dealii/dealii/pull/12848

#12847: [WIP] DataOut for pyramids: reorder at the end (proposed by peterrum) 
https://github.com/dealii/dealii/pull/12847

#12846: Rename Functions::LevelSet -> Functions::SignedDistance (proposed by 
peterrum; merged) https://github.com/dealii/dealii/pull/12846

#12845: Fix order of vertices of wedges in ReferenceCell::vertex (proposed by 
peterrum; merged) https://github.com/dealii/dealii/pull/12845

#12844: Replace uses of GeometryInfo. (proposed by bangerth; merged) 
https://github.com/dealii/dealii/pull/12844

#12843: Systematize vertex construction. (proposed by bangerth; merged) 
https://github.com/dealii/dealii/pull/12843

#12841: Simplify some code by using std::any_of(). (proposed by bangerth; 
merged) https://github.com/dealii/dealii/pull/12841

#12840: Fix MatrixCreator::create_mass/laplace_matrix for hp (proposed by 
peterrum; merged) https://github.com/dealii/dealii/pull/12840

#12839: Use some more std::arrays. (proposed by drwells; merged) 
https://github.com/dealii/dealii/pull/12839

#12838: FE_Q_iso_Q1: varying subdivisions (proposed by peterrum) 
https://github.com/dealii/dealii/pull/12838

#12837: QIterated: varying subdivisions (proposed by peterrum; merged) 
https://github.com/dealii/dealii/pull/12837

#12836: Fix determination of identity in MGTwoLevelTransfer (proposed by 
peterrum; merged) https://github.com/dealii/dealii/pull/12836

#12835: Add include to mg_transfer_global_coarsening.h (proposed by peterrum; 
merged) https://github.com/dealii/dealii/pull/12835

#12834: Add one, fix another assertion. (proposed by bangerth; merged) 
https://github.com/dealii/dealii/pull/12834

#12833: Simplify a piece of code. (proposed by bangerth; merged) 
https://github.com/dealii/dealii/pull/12833

#12832: variable checkpointing fixes for >4GB files (proposed by tjhei) 
https://github.com/dealii/dealii/pull/12832

#12830: fix checkpointing for >4GB files (proposed by tjhei; merged) 
https://github.com/dealii/dealii/pull/12830

#12829: Add assert to check the definition of CellData (proposed by peterrum; 
merged) https://github.com/dealii/dealii/pull/12829

#12822: Minor text updates to the step-49 gmsh discussion. (proposed by 
bangerth; merged) https://github.com/dealii/dealii/pull/12822

#12819: Enable PreconditionChebyshev to update a sub-range of vector entries 
(proposed by kronbichler; merged) https://github.com/dealii/dealii/pull/12819

#12816: Generalize `GridTools::rotate` for arbitrary rotation axes. (proposed 
by marcfehling; merged) https://github.com/dealii/dealii/pull/12816

#12811: Correction of comments in step-72 (proposed by gfcas; merged) 
https://github.com/dealii/dealii/pull/12811

#12787: Ignore any sections before 2411 in `GridIn::read_unv()` (proposed by 
vachan-potluri; merged) https://github.com/dealii/dealii/pull/12787

#12341: Feature p4est find partition - Find MPI ranks of point owners in 
distributed meshes (proposed by konsim83; merged) 
https://github.com/dealii/dealii/pull/12341


## And this is a list of recently opened or closed discussions:

#12854: constant_modes for PreconditionAMG (opened) 
https://github.com/dealii/dealii/issues/12854

#12842: Ordering of vertices for wedges and pyramids (opened) 
https://github.com/dealii/dealii/issues/12842

#12831: compiling issue (opened and closed) 
https://github.com/dealii/dealii/issues/12831


A list of all major changes since the last release can be found at 
https://www.dealii.org/developer/doxygen/deal.II/recent_changes.html.


Thanks for being part of the community!


Let us know about questions, problems, bugs or just share your experience by 
writing to dealii@googlegroups.com, or by opening issues or pull requests at 
https://www.github.com/dealii/dealii.
Additional information can be found at https://www.dealii.org/.

-- 
The deal.II project 

[deal.II] Preconditioner for asymmetric matrix-free style matrix

2021-10-21 Thread 殷承江
Hi all,

As described in the title, I am looking for an appropriate preconditioner 
for an asymmetric matrix implemented in matrix-free style.

I've gone though all the matrix-free tutorials and none of them talks about 
the asymmetric case. (Almost all of them use the Chebyshev smoother with 
different inner preconditioners)

However, *Chebyshev smoother* is out of the question due to its requirement 
of symmetric positive definite matrix to get the max eigen-value. 
*BlockRelaxation 
preconditioner* is also unlikely to fit, which requires the inner 
contribution to cell blocks using discontinuous finite element. Other kinds 
of relaxation preconditioner like *SOR, ILU* all request certain entries of 
the matrix, which is expensive in matrix-free framework.

So I guess the only doable choice would be the *PreconditionJacobi + GMRES 
on the inner levels and AMG + GMRES solver on the coarse grid*?

Since I am still learning this amazing library and not familiar with the 
framework. Feel free to correct me if I got anything wrong.

Regards,
Chengjiang Yin

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/0f8197d4-6b88-42d0-a2a3-303e05601a08n%40googlegroups.com.


Re: [deal.II] Dirichlet bc via apply_boundary_values function

2021-10-21 Thread Мария Бронзова
Dear Wolfgang,

I believe I answered you privately, but I cannot find my email anywhere. In 
case you didn"t receive it, I am sending such a test case here again. The 
testcase is a box porous material model with 8 cells to be able to 
interpolate BC's and to see the problematics clearlier. The confusing thing 
might be the Time function, as I am working with the frequency domain, but 
I simply used it artificially as a counter for frequencies to generate a 
series of output files. 

The problem appears in the pressure boundary condition, as I am trying to 
force (250,0) complex-valued boundary values on certain faces. After 
solving the system, the Dirichlet boundary values are not held, so (250,0) 
becomes somewhat else.

Please, let me know, if some other information is required!

Thank you a lot!

Kind regards,
Mariia



#include 
#include 
#include 

#include 
#include 
#include 
#include 
#include 

#include 
#include 

#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 

#include 
#include 

#include 
#include 
#include 
#include 
#include 

#include 
#include 
#include 
#include 

using namespace std;

namespace Poroacoustics
{
  using namespace dealii;

  template 
  SymmetricTensor<4, dim, std::complex> 
get_stress_strain_tensor(const std::complex lambda,
   const std::complex mu)
  {
  const std::complex nothing(0.0,0.0);
  SymmetricTensor<4, dim, std::complex> stress_strain_tensor;
  for (unsigned int i=0; i< dim; ++i)
  for (unsigned int j=0; j< dim; ++j)
  for (unsigned int k=0; k< dim; ++k)
  for (unsigned int l=0; l< dim; ++l)
  stress_strain_tensor[i][j][k][l]= (((i==k) && (j==l) ? mu : nothing) +
  ((i==l) && (j==k) ? mu : nothing) +
((i==j) && (k==l) ? lambda : nothing));
  return stress_strain_tensor;
  }

  template 
  inline SymmetricTensor<2, dim> get_strain(const FEValues &fe_values,
  const unsigned int shape_func,
const unsigned int q_point)
  {
  SymmetricTensor<2, dim> strain_tensor;
  for (unsigned int i=0; i< dim; ++i)
  strain_tensor[i][i] = fe_values.shape_grad_component(shape_func, q_point, 
i)[i];

  for (unsigned int i=0; i< dim; ++i)
  for (unsigned int j=i+1; j< dim; ++j)
  strain_tensor[i][j] =
  (fe_values.shape_grad_component(shape_func, q_point, i)[j] +
   fe_values.shape_grad_component(shape_func, q_point, j)[i]) /
   2;
  return strain_tensor;
  }


  template 
  class Poroelasticity
  {
  public:
Poroelasticity(const unsigned int degree);
void run();

  private:
void make_grid_and_dofs();
void assemble_system();
void solve();
void output_results() const;

const unsigned int degree;

Triangulation triangulation;
FESystem  fe;
DoFHandlerdof_handler;
DiscreteTime time;


BlockSparsityPattern  sparsity_pattern;
BlockSparseMatrix> system_matrix;
BlockVector> solution;
BlockVector> system_rhs;

static const SymmetricTensor<4, dim> stress_strain_tensor;

static constexpr unsigned int start_freq=100;
static constexpr unsigned int end_freq=101; //[Hz]

  };

  template
  class RightHandSide : public Function>
  {
  public:
  RightHandSide()
  : Function>(dim+1)
  {}

  virtual std::complex value(const Point &p,
 const unsigned int component = 0) const override;
  };

  template 
  std::complex RightHandSide::value(const Point & p,
   const unsigned int component) const
  {
  const std::complex i = {0.,0.};

  return {0.,0.};
  }


  template 
  class PressureBoundaryValues : public Function>
  {
  public:
  PressureBoundaryValues()
: Function>(dim+1)
  {}
  virtual std::complex value(const Point &p,
 const unsigned int component = 0) const override;

  };

  template 
  std::complex PressureBoundaryValues::value(const Point 
& p,
   const unsigned int component) const
  {

  const std::complex i = {0,1};
  return 250.+ 0.0*i;

  }

  template 
  class DisplBoundaryValues : public Function>
  {
  public:
  DisplBoundaryValues()
: Function>(dim+1)
  {}

  virtual std::complex value(const Point &p,
   const unsigned int component = 0) const override;

  virtual void vector_value(const Point &p,
  Vector> &  value) const override;
  };

template 
  std::complex DisplBoundaryValues::value(const Point & p,
 const unsigned int component) const
{
   //const std::complex i = {0,1};
   return {0.0,0.0};
}

template 
void DisplBoundaryValues::vector_value(const Point &p,
   Vector> & 
values) const
{

  const std::complex i = {0,1};
  values(0) = DisplBoundaryValues::value(p, 0);
  values(1) = DisplBoundaryValues::value(p, 1);
  values(2) = 0.0+i*0.0;

}



  template 
  Poroelasticity::Poroelasticity(const unsigned int degree)
: degree(degree)
  , triangulation(Triangulation::maximum_smoothing)
, fe(FE_Q(degree

[deal.II] n components with matrix-free and MPI

2021-10-21 Thread Hermes Sampedro
Dear all, 

I am having some trouble understanding how to have access to the different 
components (2 for complex values) in a matrix-free implementation and MPI. 
I need to have access to each component to perform the different operations 
to them in the *local_apply()* function.

I saw this documentation 
(https://www.dealii.org/current/doxygen/deal.II/classFEEvaluation.html, 
@Handling multi-component systems) which get the components using 
BlockVectors. 
However, I define in my class *FESystem  fe;  *and after, 
*fe(FE_Q(degree), 
2) * to have 2 components and use LinearAlgebra::distributed::Vector 

 for 
the solution, dst, scr.
In this case, is it possible to have access to the different components 
using blocks() as it is done in the example (
https://github.com/dealii/dealii/blob/master/tests/matrix_free/matrix_vector_stokes_noflux.cc,
 
line 67) using for example phi.read_dof_values(src.block(0)) ?


I would really appreciate it if I can get help on this matter.

Regards, 
H.



-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/320fbf7d-9532-4c2f-9952-46e3f501316en%40googlegroups.com.