int main(int argc, char *argv[])
{
bool plot = false;
gsCmdLine cmd(
"Tutorial on assemblying a Poisson problem.");
cmd.addSwitch("plot", "Create a ParaView visualization file with the solution", plot);
try { cmd.getValues(argc,argv); } catch (int rv) { return rv; }
gsInfo <<
"The domain is: "<< patches <<
"\n";
2);
for (gsMultiPatch<>::const_biterator
bit = patches.
bBegin(); bit != patches.
bEnd(); ++bit)
{
}
int numRefine = 2;
int numElevate = 2;
if ( numElevate > -1 )
{
int max_tmp = splinebasis.minCwiseDegree();
max_tmp += numElevate;
splinebasis.setDegree(max_tmp);
}
for (int i = 0; i < numRefine; ++i)
splinebasis.uniformRefine();
2);
opt.
setInt(
"DirichletValues" , dirichlet::l2Projection);
opt.
setInt(
"DirichletStrategy", dirichlet::elimination );
opt.
setInt(
"InterfaceStrategy", iFace ::conforming );
splinebasis.getMapper((dirichlet::strategy)opt.getInt("DirichletStrategy"),
(iFace:: strategy)opt.getInt("InterfaceStrategy"),
bcInfo, mapper, 0);
sys.reserve(10, ppde.numRhs() );
gsInfo <<
"Assembled a system (matrix and load vector) with "
gsSparseSolver<>::CGDiagonal solver( PA.
matrix() );
gsInfo <<
"Solved the system with CG solver.\n";
if (plot)
{
gsInfo<<
"Plotting in Paraview...\n";
gsWriteParaview<>(sol, "poisson2d", 1000);
gsWriteParaview<>( exact, "poisson2d_exact", 1000);
}
else
{
gsInfo <<
"Done. No output created, re-run with --plot to get a ParaView "
"file containing the solution.\n";
}
return 0;
}
The assembler class provides generic routines for volume and boundary integrals that are used for for...
Definition gsAssembler.h:266
index_t numDofs() const
Returns the number of (free) degrees of freedom.
Definition gsAssembler.h:633
void initialize(const gsPde< T > &pde, const gsStdVectorRef< gsMultiBasis< T > > &bases, const gsOptionList &opt=defaultOptions())
Intitialize function for, sets data fields using the pde, a vector of multi-basis and assembler optio...
Definition gsAssembler.h:317
const gsMatrix< T > & rhs() const
Returns the left-hand side vector(s) ( multiple right hand sides possible )
Definition gsAssembler.h:618
void finalize()
finishes the assembling of the system matrix, i.e. calls its .makeCompressed() method.
Definition gsAssembler.h:388
virtual void constructSolution(const gsMatrix< T > &solVector, gsMultiPatch< T > &result, short_t unk=0) const
Construct solution from computed solution vector for a single unknows.
Definition gsAssembler.hpp:537
void setSparseSystem(gsSparseSystem< T > &sys)
Swaps the actual sparse system with the given one.
Definition gsAssembler.h:626
void computeDirichletDofs(short_t unk=0)
Triggers computation of the Dirichlet dofs.
Definition gsAssembler.hpp:232
static gsOptionList defaultOptions()
Returns the list of default options for assembly.
Definition gsAssembler.hpp:30
void push()
Iterates over all elements of the domain and applies the ElementVisitor.
Definition gsAssembler.h:428
const gsMultiPatch< T > & patches() const
Return the multipatch.
Definition gsAssembler.h:601
const gsSparseMatrix< T > & matrix() const
Returns the left-hand global matrix.
Definition gsAssembler.h:614
Class containing a set of boundary conditions.
Definition gsBoundaryConditions.h:342
void addCondition(int p, boxSide s, condition_type::type t, gsFunctionSet< T > *f, short_t unknown=0, bool parametric=false, int comp=-1)
Adds another boundary condition.
Definition gsBoundaryConditions.h:650
const_biterator bBegin() const
Definition gsBoxTopology.h:139
const_biterator bEnd() const
Definition gsBoxTopology.h:144
Class for command-line argument parsing.
Definition gsCmdLine.h:57
Maintains a mapping from patch-local dofs to global dof indices and allows the elimination of individ...
Definition gsDofMapper.h:69
std::ostream & print(std::ostream &os=gsInfo) const
Print summary.
Definition gsDofMapper.cpp:346
A scalar of vector field defined on a m_parametric geometry.
Definition gsField.h:55
static void open(const std::string &fn)
Opens the file fn using the preferred application of the OS.
Definition gsFileManager.cpp:688
Class defining a multivariate (real or vector) function given by a string mathematical expression.
Definition gsFunctionExpr.h:52
A matrix with arbitrary coefficient type and fixed or dynamic size.
Definition gsMatrix.h:41
Holds a set of patch-wise bases and their topology information.
Definition gsMultiBasis.h:37
Container class for a set of geometry patches and their topology, that is, the interface connections ...
Definition gsMultiPatch.h:100
Class which holds a list of parameters/options, and provides easy access to them.
Definition gsOptionList.h:33
void setReal(const std::string &label, const Real &value)
Sets an existing option label to be equal to value.
Definition gsOptionList.cpp:166
void setInt(const std::string &label, const index_t &value)
Sets an existing option label to be equal to value.
Definition gsOptionList.cpp:158
A Poisson PDE.
Definition gsPoissonPde.h:35
A sparse linear system indexed by sets of degrees of freedom.
Definition gsSparseSystem.h:30
Implementation of a Neumann BC for elliptic assemblers.
Definition gsVisitorNeumann.h:34
Visitor for the Poisson equation.
Definition gsVisitorPoisson.h:35
Main header to be included by clients using the G+Smo library.
Provides generic assembler routines.
#define gsInfo
Definition gsDebug.h:43
Neumann conditions visitor for elliptic problems.
Weak (Nitsche-type) imposition of the Dirichlet boundary conditions.
Poisson equation element visitor.
The G+Smo namespace, containing all definitions for the library.
@ dirichlet
Dirichlet type.
Definition gsBoundaryConditions.h:31
static gsMultiPatch< T > BSplineSquareGrid(int n, int m, T const &r=1, T const &lx=0, T const &ly=0)
Definition gsNurbsCreator.hpp:691