Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

DAE initialization note #604

Merged
merged 6 commits into from
Feb 28, 2023
Merged
Show file tree
Hide file tree
Changes from 2 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions docs/make.jl
Original file line number Diff line number Diff line change
Expand Up @@ -52,7 +52,7 @@ GENERATEDEXAMPLES = [joinpath("examples", f) for f in (
"Manual" => [
"manual/degrees_of_freedom.md",
"manual/assembly.md",
"manual/boundary_conditions.md",
"manual/initial_and_boundary_conditions.md",
fredrikekre marked this conversation as resolved.
Show resolved Hide resolved
"manual/constraints.md",
"manual/grid.md",
"manual/export.md"
Expand All @@ -64,7 +64,7 @@ GENERATEDEXAMPLES = [joinpath("examples", f) for f in (
"reference/fevalues.md",
"reference/dofhandler.md",
"reference/assembly.md",
"reference/boundary_conditions.md",
"reference/initial_and_boundary_conditions.md",
fredrikekre marked this conversation as resolved.
Show resolved Hide resolved
"reference/grid.md",
"reference/export.md",
"reference/utils.md",
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
DocTestSetup = :(using Ferrite)
```

# Boundary Conditions
# Initial and Boundary Conditions

Every PDE is accompanied with boundary conditions. There are different types of boundary
conditions, and they need to be handled in different ways. Below we discuss how to handle
Expand Down Expand Up @@ -256,7 +256,7 @@ pdbc = PeriodicDirichlet(
)
```

# Initial Conditions
## Initial Conditions

When solving time-dependent problems, initial conditions, different from zero, may be required.
For finite element formulations of ODE-type,
Expand All @@ -272,4 +272,32 @@ u = zeros(ndofs(dh))
apply_analytical!(u, dh, :p, x -> ρ * g * x[2])
```

!!! note "Initialization of Differential-Algebraic Systems"
Differential-algebraic systems of equations (DAEs) need extra care during initialization, which is
currently not automatized by the function `apply_analytical!`.
Skipping the formal definition, you can identify DAEs for example by
* Multiple time derivatives in the equation
* The equations cannot be formulated such that the variable with the time-derivative
is the only term on the right hand side
* Missing time derivatives
To give concrete example of a "DAE" (formally a PDAE, which becomes a DAE by spatial discretization), consider following
abstract problem on some domain with ``u(t,\boldsymbol{x})`` and ``s(t,\boldsymbol{x})`` being the unknown variables:
```math
\partial_t u(t,\boldsymbol{x})) = f(u(t,\boldsymbol{x}), \nabla u(t,\boldsymbol{x}), s(t,\boldsymbol{x})) \\
0 = g(u(t,\boldsymbol{x}), \nabla u(t,\boldsymbol{x}), s(t,\boldsymbol{x}))
```
Here ``s`` might be some internal state variable, e.g. describing hardening, failure or some electrical
state and ``u`` could be a displacement field or some potential field. This form is quite common in mechanics,
for example in quasi-static analysis of plastic materials.
For DAEs need to assert that the hidden constraints of the problem are not violated by the choice
of initial conditions.
What the exact hidden constraints are have to be worked out on the specific DAE.
For the example above it is sometimes possible to just solve the second equation for ``s``, i.e.
```math
0 = g(u(t,\boldsymbol{x}), \nabla u(t,\boldsymbol{x}), s(t,\boldsymbol{x}))
```
given some fixed choice of ``u``.
We refer to the paper ["Consistent Initial Condition Calculation for Differential-Algebraic Systems"
by Brown et al.](dx.doi.org/10.1137/S1064827595289996) for more details on this matter and some possible solutions.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for providing such a detailed description Dennis, it is very interesting!
But this does indeed seem very difficult to explain in simple terms for the general case to not confuse (new) users.
Thinking more about this, could something like the following, even if over-simplistic, work?

Suggested change
!!! note "Initialization of Differential-Algebraic Systems"
Differential-algebraic systems of equations (DAEs) need extra care during initialization, which is
currently not automatized by the function `apply_analytical!`.
Skipping the formal definition, you can identify DAEs for example by
* Multiple time derivatives in the equation
* The equations cannot be formulated such that the variable with the time-derivative
is the only term on the right hand side
* Missing time derivatives
To give concrete example of a "DAE" (formally a PDAE, which becomes a DAE by spatial discretization), consider following
abstract problem on some domain with ``u(t,\boldsymbol{x})`` and ``s(t,\boldsymbol{x})`` being the unknown variables:
```math
\partial_t u(t,\boldsymbol{x})) = f(u(t,\boldsymbol{x}), \nabla u(t,\boldsymbol{x}), s(t,\boldsymbol{x})) \\
0 = g(u(t,\boldsymbol{x}), \nabla u(t,\boldsymbol{x}), s(t,\boldsymbol{x}))
```
Here ``s`` might be some internal state variable, e.g. describing hardening, failure or some electrical
state and ``u`` could be a displacement field or some potential field. This form is quite common in mechanics,
for example in quasi-static analysis of plastic materials.
For DAEs need to assert that the hidden constraints of the problem are not violated by the choice
of initial conditions.
What the exact hidden constraints are have to be worked out on the specific DAE.
For the example above it is sometimes possible to just solve the second equation for ``s``, i.e.
```math
0 = g(u(t,\boldsymbol{x}), \nabla u(t,\boldsymbol{x}), s(t,\boldsymbol{x}))
```
given some fixed choice of ``u``.
We refer to the paper ["Consistent Initial Condition Calculation for Differential-Algebraic Systems"
by Brown et al.](dx.doi.org/10.1137/S1064827595289996) for more details on this matter and some possible solutions.
!!! note "Initialization of Differential-Algebraic Systems"
Note that `apply_analytical!` does not ensure consistency between the applied field
and the specific problem being solved. As a simple example, consider a mechanical problem:
`apply_analytical!` doesn't ensure that equilibrium is fulfilled in the initial state.
In general for advanced, so-called "Differential-Algebraic Systems",
consistent initialization is a research field on its own, see e.g.
["Consistent Initial Condition Calculation for Differential-Algebraic Systems" by Brown et al.](dx.doi.org/10.1137/S1064827595289996)

or do you think this misses the point since in many cases a non-equilibrium initial condition is perfectly fine (and we even demo that in the transient heat flow:))

I am, however, sad not to have all such knowledge/insights available. From that perspective, it might have been nice to have a separate webpage containing "Ferrite: Theory and examples" (or something) to supplement the documentation, but that's a different discussion

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for looking into this and the suggestion. I am happy and open for ideas on how to improve/simplify this part. However, the point here is not about equilibrium, but about consistency of the solution. This is especially problematic in differential-algebraic problems.

Think about it this way: The algebraic equations act as constraints on the time evolution of your problem, where the constraints also have some "implicit" or "hidden evolution". However, the constraint might not be directly visible in the equation and hence messing it up can be very easy (for example if a zero initialization of some internal variable could be inconsistent with your evolving field variable).

Don't you think the (theory) manual we have is sufficient in combination with the examples?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, I agree that equilibrium is not the point. I meant it more as a warning that one can not supply whatever initial condition but have to ensure that it is consistent with the problem to be solved.

(The additional site was more a long term idea that would be useful for e.g. teaching and reference since you, we and maybe more use Ferrite in teaching, for which I think the current state is not meant to be, nor detailed enough, as course material)

See also [Time Dependent Problems](@ref) for one example.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should we put this above the note?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Quick question here, why would you move it above the note?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It is easy to miss that reference when reading the docs, but if switched, the note is still hard to miss. And it follows (in my mind) more logically as DAEs are not discussed in that example.

Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@ DocTestSetup = :(using Ferrite)
# Boundary Conditions

```@index
Pages = ["boundary_conditions.md"]
Pages = ["initial_and_boundary_conditions.md"]
fredrikekre marked this conversation as resolved.
Show resolved Hide resolved
```

```@docs
Expand Down