Skip to content
GitLab
Projects
Groups
Snippets
/
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
amdis
amdis-core
Commits
d9715c4a
Commit
d9715c4a
authored
Jan 31, 2020
by
Müller, Felix
Browse files
Add overload for LagrangeBasis with single node
parent
eda87e6c
Changes
18
Hide whitespace changes
Inline
Side-by-side
examples/CMakeLists.txt
View file @
d9715c4a
...
...
@@ -35,4 +35,4 @@ endif ()
if
(
dune-foamgrid_FOUND
)
add_amdis_executable
(
NAME surface.2d SOURCES surface.cc DIM 2 DOW 3 ALBERTA_GRID
)
add_dependencies
(
examples surface.2d
)
endif
()
\ No newline at end of file
endif
()
examples/boundary.cc
View file @
d9715c4a
...
...
@@ -11,10 +11,8 @@
#include
<amdis/Integrate.hpp>
#include
<amdis/ProblemStat.hpp>
#include
<amdis/LocalOperators.hpp>
#include
<amdis/common/Literals.hpp>
using
namespace
AMDiS
;
using
namespace
Dune
::
Indices
;
// 1 component with polynomial degree 2
using
Param
=
YaspGridBasis
<
GRIDDIM
,
2
>
;
...
...
@@ -28,7 +26,7 @@ void run(SetBoundary setBoundary)
setBoundary
(
prob
.
boundaryManager
());
auto
opL
=
makeOperator
(
tag
::
gradtest_gradtrial
{},
1.0
);
prob
.
addMatrixOperator
(
opL
,
0
,
0
);
prob
.
addMatrixOperator
(
opL
);
auto
f
=
[](
auto
const
&
x
)
{
double
r2
=
dot
(
x
,
x
);
...
...
@@ -36,18 +34,18 @@ void run(SetBoundary setBoundary)
return
-
(
400.0
*
r2
-
20.0
*
x
.
size
())
*
ux
;
};
auto
opForce
=
makeOperator
(
tag
::
test
{},
f
,
6
);
prob
.
addVectorOperator
(
opForce
,
0
);
prob
.
addVectorOperator
(
opForce
);
// set boundary condition
auto
g
=
[](
auto
const
&
x
){
return
std
::
exp
(
-
10.0
*
dot
(
x
,
x
));
};
prob
.
addDirichletBC
(
BoundaryType
{
1
},
0
,
0
,
g
);
prob
.
addDirichletBC
(
BoundaryType
{
1
},
g
);
AdaptInfo
adaptInfo
(
"adapt"
);
prob
.
assemble
(
adaptInfo
);
prob
.
solve
(
adaptInfo
);
double
errorL2
=
integrate
(
sqr
(
g
-
prob
.
solution
(
0
)),
prob
.
gridView
(),
6
);
double
errorL2
=
integrate
(
sqr
(
g
-
prob
.
solution
()),
prob
.
gridView
(),
6
);
msg
(
"error_L2 = {}"
,
errorL2
);
}
...
...
@@ -69,17 +67,17 @@ void run_periodic()
prob
.
boundaryManager
().
setBoxBoundary
({
-
1
,
-
1
,
1
,
1
});
auto
opL
=
makeOperator
(
tag
::
gradtest_gradtrial
{},
1.0
);
prob
.
addMatrixOperator
(
opL
,
0
,
0
);
prob
.
addMatrixOperator
(
opL
);
auto
opForce
=
makeOperator
(
tag
::
test
{},
1.0
,
6
);
prob
.
addVectorOperator
(
opForce
,
0
);
prob
.
addVectorOperator
(
opForce
);
// set boundary condition
auto
g
=
[](
auto
const
&
x
)
{
return
std
::
sin
(
0.5
*
M_PI
*
x
[
1
])
*
std
::
sin
(
2
*
M_PI
*
(
0.25
+
x
[
0
]))
+
std
::
cos
(
0.5
*
M_PI
*
x
[
1
])
*
std
::
sin
(
2
*
M_PI
*
(
-
0.25
+
x
[
0
]));
};
prob
.
addDirichletBC
(
BoundaryType
{
1
},
0
,
0
,
g
);
prob
.
addDirichletBC
(
BoundaryType
{
1
},
g
);
typename
ElliptProblem
::
WorldMatrix
A
{{
1.0
,
0.0
},
{
0.0
,
1.0
}};
typename
ElliptProblem
::
WorldVector
b
{
1.0
,
0.0
};
...
...
examples/convection_diffusion.cc
View file @
d9715c4a
...
...
@@ -6,7 +6,6 @@
#include
<amdis/AMDiS.hpp>
#include
<amdis/ProblemStat.hpp>
#include
<amdis/localoperators/ConvectionDiffusionOperator.hpp>
#include
<amdis/common/Literals.hpp>
using
namespace
AMDiS
;
...
...
@@ -19,20 +18,18 @@ int main(int argc, char** argv)
{
Environment
env
(
argc
,
argv
);
using
namespace
Dune
::
Indices
;
ElliptProblem
prob
(
"ellipt"
);
prob
.
initialize
(
INIT_ALL
);
// -div(A*grad(u)) + div(b*u) + c*u = f
auto
opCD
=
convectionDiffusion
(
/*A=*/
1.0
,
/*b=*/
0.0
,
/*c=*/
1.0
,
/*f=*/
1.0
);
prob
.
addMatrixOperator
(
opCD
,
0
,
0
);
prob
.
addVectorOperator
(
opCD
,
0
);
prob
.
addMatrixOperator
(
opCD
);
prob
.
addVectorOperator
(
opCD
);
// set boundary condition
auto
predicate
=
[](
auto
const
&
x
){
return
x
[
0
]
<
1.e-8
||
x
[
1
]
<
1.e-8
;
};
// define boundary
auto
dbcValues
=
[](
auto
const
&
x
){
return
0.0
;
};
// set value
prob
.
addDirichletBC
(
predicate
,
0
,
0
,
dbcValues
);
prob
.
addDirichletBC
(
predicate
,
dbcValues
);
AdaptInfo
adaptInfo
(
"adapt"
);
prob
.
assemble
(
adaptInfo
);
...
...
examples/ellipt.cc
View file @
d9715c4a
...
...
@@ -8,7 +8,6 @@
using
Grid
=
Dune
::
YaspGrid
<
GRIDDIM
>
;
using
namespace
AMDiS
;
using
namespace
Dune
::
Indices
;
int
main
(
int
argc
,
char
**
argv
)
{
...
...
@@ -34,13 +33,13 @@ int main(int argc, char** argv)
prob
.
initialize
(
INIT_ALL
);
auto
opL
=
makeOperator
(
tag
::
gradtest_gradtrial
{},
1.0
);
prob
.
addMatrixOperator
(
opL
,
0
,
0
);
prob
.
addMatrixOperator
(
opL
);
auto
opForce
=
makeOperator
(
tag
::
test
{},
f
,
6
);
prob
.
addVectorOperator
(
opForce
,
0
);
prob
.
addVectorOperator
(
opForce
);
// set boundary condition
prob
.
addDirichletBC
(
1
,
0
,
0
,
g
);
prob
.
addDirichletBC
(
1
,
g
);
AdaptInfo
adaptInfo
(
"adapt"
);
...
...
@@ -65,9 +64,9 @@ int main(int argc, char** argv)
prob
.
assemble
(
adaptInfo
);
prob
.
solve
(
adaptInfo
);
double
errorL2
=
integrate
(
sqr
(
g
-
prob
.
solution
(
0
)),
prob
.
gridView
(),
6
);
double
errorL2
=
integrate
(
sqr
(
g
-
prob
.
solution
()),
prob
.
gridView
(),
6
);
errL2
.
push_back
(
std
::
sqrt
(
errorL2
));
double
errorH1
=
errorL2
+
integrate
(
unary_dot
(
grad_g
-
gradientAtQP
(
prob
.
solution
(
0
))),
prob
.
gridView
(),
6
);
double
errorH1
=
errorL2
+
integrate
(
unary_dot
(
grad_g
-
gradientAtQP
(
prob
.
solution
())),
prob
.
gridView
(),
6
);
errH1
.
push_back
(
std
::
sqrt
(
errorH1
));
}
...
...
examples/heat.cc
View file @
d9715c4a
...
...
@@ -9,7 +9,6 @@
#include
<amdis/ProblemInstat.hpp>
#include
<amdis/ProblemStat.hpp>
#include
<amdis/GridFunctions.hpp>
#include
<amdis/common/Literals.hpp>
using
namespace
AMDiS
;
...
...
@@ -34,23 +33,23 @@ int main(int argc, char** argv)
auto
invTau
=
std
::
ref
(
probInstat
.
invTau
());
auto
opTimeLhs
=
makeOperator
(
tag
::
test_trial
{},
invTau
);
prob
.
addMatrixOperator
(
opTimeLhs
,
0
,
0
);
prob
.
addMatrixOperator
(
opTimeLhs
);
auto
opL
=
makeOperator
(
tag
::
gradtest_gradtrial
{},
1.0
);
prob
.
addMatrixOperator
(
opL
,
0
,
0
);
prob
.
addMatrixOperator
(
opL
);
auto
opTimeRhs
=
makeOperator
(
tag
::
test
{},
invokeAtQP
([
invTau
](
double
u
)
{
return
u
*
invTau
.
get
();
},
prob
.
solution
(
0
)),
2
);
prob
.
addVectorOperator
(
opTimeRhs
,
0
);
invokeAtQP
([
invTau
](
double
u
)
{
return
u
*
invTau
.
get
();
},
prob
.
solution
()),
2
);
prob
.
addVectorOperator
(
opTimeRhs
);
auto
opForce
=
makeOperator
(
tag
::
test
{},
[](
auto
const
&
x
)
{
return
-
1.0
;
},
0
);
prob
.
addVectorOperator
(
opForce
,
0
);
prob
.
addVectorOperator
(
opForce
);
// set boundary condition
auto
predicate
=
[](
auto
const
&
p
){
return
p
[
0
]
<
1.e-8
||
p
[
1
]
<
1.e-8
;
};
auto
dbcValues
=
[](
auto
const
&
p
){
return
0.0
;
};
prob
.
addDirichletBC
(
predicate
,
0
,
0
,
dbcValues
);
prob
.
addDirichletBC
(
predicate
,
dbcValues
);
AdaptInstationary
adapt
(
"adapt"
,
prob
,
adaptInfo
,
probInstat
,
adaptInfo
);
adapt
.
adapt
();
...
...
examples/init/boundary.dat.2d
View file @
d9715c4a
...
...
@@ -8,7 +8,7 @@ ellipt->solver->relative tolerance: 1.e-10
ellipt->solver->info: 1
ellipt->solver->precon: ilu
ellipt->output
[0]
->format: vtk
ellipt->output
[0]
->filename: boundary.2d
ellipt->output
[0]
->name: u
ellipt->output
[0]
->output directory: ./output
ellipt->output->format: vtk
ellipt->output->filename: boundary.2d
ellipt->output->name: u
ellipt->output->output directory: ./output
examples/init/cahn_hilliard.dat.2d
View file @
d9715c4a
...
...
@@ -16,7 +16,7 @@ ch->solver: direct
ch->solver->max iteration: 1000
ch->solver->relative tolerance: 1e-6
ch->solver->break if tolerance not reached: 1
ch->solver->info:
2
ch->solver->info:
0
ch->output[0]->format: vtk
ch->output[0]->filename: phi.2d
...
...
@@ -34,4 +34,4 @@ ch->output[1]->animation: 1
adapt->timestep: 0.001
adapt->start time: 0.0
adapt->end time:
1.0
adapt->end time:
0.5
examples/init/ellipt.dat.2d
View file @
d9715c4a
...
...
@@ -13,7 +13,7 @@ ellipt->solver->max iteration: 10000
ellipt->solver->relative tolerance: 1.e-10
ellipt->solver->precon: ilu
ellipt->output
[0]
->format: vtk backup
ellipt->output
[0]
->filename: ellipt.2d
ellipt->output
[0]
->name: u
ellipt->output
[0]
->output directory: ./output
ellipt->output->format: vtk backup
ellipt->output->filename: ellipt.2d
ellipt->output->name: u
ellipt->output->output directory: ./output
examples/init/ellipt.dat.3d
View file @
d9715c4a
...
...
@@ -13,7 +13,7 @@ ellipt->solver->max iteration: 10000
ellipt->solver->relative tolerance: 1.e-10
ellipt->solver->precon: ilu
ellipt->output
[0]
->format: vtk
ellipt->output
[0]
->filename: ellipt.3d
ellipt->output
[0]
->name: u
ellipt->output
[0]
->output directory: ./output
ellipt->output->format: vtk
ellipt->output->filename: ellipt.3d
ellipt->output->name: u
ellipt->output->output directory: ./output
examples/init/heat.dat.2d
View file @
d9715c4a
...
...
@@ -12,12 +12,12 @@ heat->solver->absolute tolerance: 1e-6
heat->solver->info: 1
heat->solver->precon: ilu
heat->output
[0]
->format: vtk backup
heat->output
[0]
->name: u
heat->output
[0]
->filename: heat.2d
heat->output
[0]
->output directory: output
heat->output
[0]
->mode: 1
heat->output
[0]
->animation: 1
heat->output->format: vtk backup
heat->output->name: u
heat->output->filename: heat.2d
heat->output->output directory: output
heat->output->mode: 1
heat->output->animation: 1
adapt->timestep: 0.1
adapt->start time: 0.0
...
...
examples/init/surface.dat.2d
View file @
d9715c4a
...
...
@@ -9,7 +9,7 @@ surface->solver->relative tolerance: 1.e-6
surface->solver->info: 1
surface->solver->precon: ilu
surface->output
[0]
->format: vtk
surface->output
[0]
->filename: surface.2d
surface->output
[0]
->name: u
surface->output
[0]
->output directory: ./output
surface->output->format: vtk
surface->output->filename: surface.2d
surface->output->name: u
surface->output->output directory: ./output
examples/neumann.cc
View file @
d9715c4a
...
...
@@ -15,10 +15,8 @@
#include
<amdis/Integrate.hpp>
#include
<amdis/ProblemStat.hpp>
#include
<amdis/LocalOperators.hpp>
#include
<amdis/common/Literals.hpp>
using
namespace
AMDiS
;
using
namespace
Dune
::
Indices
;
template
<
class
Grid
>
void
run
(
Grid
&
grid
)
...
...
@@ -28,18 +26,18 @@ void run(Grid& grid)
using
Traits
=
LagrangeBasis
<
Grid
,
2
>
;
ProblemStat
<
Traits
>
prob
(
"ellipt"
,
grid
);
prob
.
initialize
(
INIT_ALL
);
prob
.
boundaryManager
()
.
setBoxBoundary
({
1
,
2
,
2
,
1
});
prob
.
boundaryManager
()
->
setBoxBoundary
({
1
,
2
,
2
,
1
});
auto
opL
=
makeOperator
(
tag
::
gradtest_gradtrial
{},
1.0
);
prob
.
addMatrixOperator
(
opL
,
0
,
0
);
prob
.
addMatrixOperator
(
opL
);
// set dirichlet boundary condition
auto
g
=
[](
auto
const
&
x
)
{
return
0.0
;
};
prob
.
addDirichletBC
(
BoundaryType
{
1
},
0
,
0
,
g
);
prob
.
addDirichletBC
(
BoundaryType
{
1
},
g
);
// set neumann boundary condition
auto
opNeumann
=
makeOperator
(
tag
::
test
{},
1.0
);
prob
.
addVectorOperator
(
BoundaryType
{
2
},
opNeumann
,
0
);
prob
.
addVectorOperator
(
BoundaryType
{
2
},
opNeumann
);
AdaptInfo
adaptInfo
(
"adapt"
);
...
...
examples/periodic.cc
View file @
d9715c4a
...
...
@@ -11,7 +11,6 @@
#include
<amdis/Integrate.hpp>
#include
<amdis/ProblemStat.hpp>
#include
<amdis/LocalOperators.hpp>
#include
<amdis/common/Literals.hpp>
using
namespace
AMDiS
;
...
...
examples/surface.cc
View file @
d9715c4a
...
...
@@ -3,7 +3,6 @@
#include
<amdis/Integrate.hpp>
#include
<amdis/LocalOperators.hpp>
#include
<amdis/ProblemStat.hpp>
#include
<amdis/common/Literals.hpp>
#include
<dune/grid/albertagrid/albertareader.hh>
#include
<dune/foamgrid/foamgrid.hh>
...
...
@@ -12,7 +11,6 @@
#include
"SphereMapping.hpp"
using
namespace
AMDiS
;
using
namespace
Dune
::
Indices
;
using
Grid
=
Dune
::
FoamGrid
<
2
,
3
>
;
using
Basis
=
LagrangeBasis
<
Grid
,
1
>
;
...
...
@@ -44,14 +42,14 @@ int main(int argc, char** argv)
prob
.
initialize
(
INIT_ALL
);
auto
opL
=
makeOperator
(
tag
::
gradtest_gradtrial
{},
1.0
);
prob
.
addMatrixOperator
(
opL
,
0
,
0
);
prob
.
addMatrixOperator
(
opL
);
double
kappa
=
Parameters
::
get
<
double
>
(
"surface->kappa"
).
value_or
(
1.0
);
auto
opM
=
makeOperator
(
tag
::
test_trial
{},
-
kappa
*
kappa
);
prob
.
addMatrixOperator
(
opM
,
0
,
0
);
prob
.
addMatrixOperator
(
opM
);
auto
opForce
=
makeOperator
(
tag
::
test
{},
X
(
0
)
+
X
(
1
)
+
X
(
2
));
prob
.
addVectorOperator
(
opForce
,
0
);
prob
.
addVectorOperator
(
opForce
);
AdaptInfo
adaptInfo
(
"adapt"
);
prob
.
assemble
(
adaptInfo
);
...
...
src/amdis/ProblemStat.hpp
View file @
d9715c4a
...
...
@@ -266,6 +266,12 @@ namespace AMDiS
RowTreePath
row
,
ColTreePath
col
,
Values
const
&
values
);
template
<
class
Identifier
,
class
Values
>
void
addDirichletBC
(
Identifier
&&
id
,
Values
&&
values
)
{
addDirichletBC
(
FWD
(
id
),
RootTreePath
{},
RootTreePath
{},
FWD
(
values
));
}
/// Add a periodic boundary conditions to the system, by specifying a face transformation
/// y = A*x + b of coordinates. We assume, that A is orthonormal.
void
addPeriodicBC
(
BoundaryType
id
,
WorldMatrix
const
&
A
,
WorldVector
const
&
b
);
...
...
src/amdis/ProblemStatTraits.hpp
View file @
d9715c4a
...
...
@@ -18,6 +18,17 @@ namespace AMDiS
template
<
bool
same
,
int
...
degs
>
struct
LagrangePreBasisCreatorImpl
;
// specialization for single node basis
template
<
int
deg
>
struct
LagrangePreBasisCreatorImpl
<
true
,
deg
>
{
static
auto
create
()
{
using
namespace
Dune
::
Functions
::
BasisFactory
;
return
lagrange
<
deg
>
();
}
};
// specialization for composite basis
template
<
int
...
degs
>
struct
LagrangePreBasisCreatorImpl
<
false
,
degs
...
>
...
...
test/ExpressionsTest.cpp
View file @
d9715c4a
...
...
@@ -7,7 +7,6 @@
#include
<amdis/Integrate.hpp>
#include
<amdis/LocalOperators.hpp>
#include
<amdis/ProblemStat.hpp>
#include
<amdis/common/Literals.hpp>
#include
<amdis/common/FieldMatVec.hpp>
...
...
@@ -20,12 +19,10 @@ int main(int argc, char** argv)
{
Environment
env
(
argc
,
argv
);
using
namespace
Dune
::
Indices
;
ElliptProblem
prob
(
"ellipt"
);
prob
.
initialize
(
INIT_ALL
);
auto
u
=
prob
.
solution
(
_0
);
auto
u
=
prob
.
solution
();
// eval a functor at global coordinates (at quadrature points)
auto
expr1
=
evalAtQP
([](
Dune
::
FieldVector
<
double
,
2
>
const
&
x
)
{
return
x
[
0
]
+
x
[
1
];
});
...
...
@@ -42,7 +39,7 @@ int main(int argc, char** argv)
auto
expr8
=
X
(
0
);
// Evaluation of the DOFVector (component) at quadrature points
auto
expr9
=
prob
.
solution
(
_0
);
auto
expr9
=
prob
.
solution
();
// ---------------------------------
...
...
test/IntegrateTest.cpp
View file @
d9715c4a
...
...
@@ -23,7 +23,7 @@ int main(int argc, char** argv)
ElliptProblem
prob
(
"ellipt"
);
prob
.
initialize
(
INIT_ALL
);
auto
u
=
prob
.
solution
(
0
);
auto
u
=
prob
.
solution
();
auto
gv
=
u
.
basis
().
gridView
();
u
<<
1.0
;
...
...
Write
Preview
Supports
Markdown
0%
Try again
or
attach a new file
.
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment