MFEM  v4.6.0 Finite element discretization library
ex37p.cpp File Reference

Go to the source code of this file.

## Functions

double proj (ParGridFunction &psi, double target_volume, double tol=1e-12, int max_its=10)
Bregman projection of ρ = sigmoid(ψ) onto the subspace ∫_Ω ρ dx = θ vol(Ω) as follows: More...

int main (int argc, char *argv[])

## ◆ main()

 int main ( int argc, char * argv[] )

## ALGORITHM PREAMBLE

The Lagrangian for this problem is

    L(u,ρ,ρ̃,w,w̃) = (f,u) - (r(ρ̃) C ε(u),ε(w)) + (f,w)
- (ϵ² ∇ρ̃,∇w̃) - (ρ̃,w̃) + (ρ,w̃)


where

r(ρ̃) = ρ₀ + ρ̃³ (1 - ρ₀) (SIMP rule)

ε(u) = (∇u + ∇uᵀ)/2 (symmetric gradient)

C e = λtr(e)I + 2μe (isotropic material)

NOTE: The Lame parameters can be computed from Young's modulus E and Poisson's ratio ν as follows:

 λ = E ν/((1+ν)(1-2ν)),      μ = E/(2(1+ν))


Discretization choices:

u ∈ V ⊂ (H¹)ᵈ (order p) ψ ∈ L² (order p - 1), ρ = sigmoid(ψ) ρ̃ ∈ H¹ (order p) w ∈ V (order p) w̃ ∈ H¹ (order p)

## ALGORITHM

Update ρ with projected mirror descent via the following algorithm.

1. Initialize ψ = inv_sigmoid(vol_fraction) so that ∫ sigmoid(ψ) = θ vol(Ω)

While not converged:

1. Solve filter equation ∂_w̃ L = 0; i.e.,

(ϵ² ∇ ρ̃, ∇ v ) + (ρ̃,v) = (ρ,v) ∀ v ∈ H¹.

2. Solve primal problem ∂_w L = 0; i.e.,

(λ r(ρ̃) ∇⋅u, ∇⋅v) + (2 μ r(ρ̃) ε(u), ε(v)) = (f,v) ∀ v ∈ V.

NB. The dual problem ∂_u L = 0 is the negative of the primal problem due to symmetry.

1. Solve for filtered gradient ∂_ρ̃ L = 0; i.e.,

(ϵ² ∇ w̃ , ∇ v ) + (w̃ ,v) = (-r'(ρ̃) ( λ |∇⋅u|² + 2 μ |ε(u)|²),v) ∀ v ∈ H¹.

2. Project the gradient onto the discrete latent space; i.e., solve
             (G,v) = (w̃,v)   ∀ v ∈ L².

3. Bregman proximal gradient update; i.e.,
                ψ ← ψ - αG + c,


where α > 0 is a step size parameter and c ∈ R is a constant ensuring

            ∫_Ω sigmoid(ψ - αG + c) dx = θ vol(Ω).


end

Definition at line 182 of file ex37p.cpp.

## ◆ proj()

 double proj ( ParGridFunction & psi, double target_volume, double tol = 1e-12, int max_its = 10 )

Bregman projection of ρ = sigmoid(ψ) onto the subspace ∫_Ω ρ dx = θ vol(Ω) as follows:

1. Compute the root of the R → R function f(c) = ∫_Ω sigmoid(ψ + c) dx - θ vol(Ω)
2. Set ψ ← ψ + c.
Parameters
 psi a GridFunction to be updated target_volume θ vol(Ω) tol Newton iteration tolerance max_its Newton maximum iteration number
Returns
double Final volume, ∫_Ω sigmoid(ψ)

Definition at line 71 of file ex37p.cpp.