Library FEM.FE
This file is part of the Coq Numerical Analysis library
Copyright (C) Boldo, Clément, Martin, Mayero, Mouhcine
This library is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public
License as published by the Free Software Foundation; either
version 3 of the License, or (at your option) any later version.
This library is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
COPYING file for more details.
[Ern]
Alexandre Ern, Éléments finis, Aide-mémoire, Dunod/L'Usine Nouvelle, 2005,
https://www.dunod.com/sciences-techniques/aide-memoire-elements-finis.
[ErnGuermond]
Alexandre Ern and Jean-Luc Guermond,
Finite Elements I. Approximation and Interpolation,
Texts in Applied Mathematics, vol. 72, Springer Cham, 2021,
https://doi.org/10.1007/978-3-030-56341-7,
https://inria.hal.science/hal-03226049 (with different page numbers).
[RR9557v1]
François Clément and Vincent Martin,
Finite element method. Detailed proofs to be formalized in Coq,
RR-9557, first version, Inria Paris, 2024,
https://inria.hal.science/hal-04713897v1.
Bibliography
https://inria.hal.science/hal-03226049 (with different page numbers).
From Requisite Require Import stdlib_wR ssr_wMC.
From Coquelicot Require Import Hierarchy.
Close Scope R_scope.
From Algebra Require Import Algebra_wDep.
From FEM Require Import geometry.
Local Open Scope R_scope.
Section FE_Facts1.
Inductive shape_type := Simplex | Quad.
Definition nvtx_of_shape (d : nat) (shape : shape_type) : nat :=
match shape with
| Simplex ⇒ d.+1
| Quad ⇒ Nat.pow 2 d
end.
Lemma nos_eq :
∀ {d1 d2 s1 s2},
d1 = d2 → s1 = s2 → nvtx_of_shape d1 s1 = nvtx_of_shape d2 s2.
Proof. move=>>; apply f_equal2. Qed.
Record FE : Type := mk_FE {
d : nat ;
ndof : nat ;
d_pos : (0 < d)%coq_nat ;
ndof_pos : (0 < ndof)%coq_nat ;
shape : shape_type ;
nvtx : nat := nvtx_of_shape d shape ;
vtx : '('R^d)^nvtx ;
K_geom : 'R^d → Prop := convex_envelop vtx ;
P_approx : FRd d → Prop ;
P_approx_has_dim : has_dim P_approx ndof ;
Sigma : '(FRd d → R)^ndof ;
Sigma_lm : ∀ i, lin_map (Sigma i) ;
unisolvence_inj : KerS0 P_approx (gather Sigma) ;
}.
Lemma unisolvence : ∀ (f:FE),
bijS (P_approx f) fullset (gather (Sigma f)).
Proof.
intros f.
apply lmS_bijS_KerS0 with (ndof f) (ndof f); try easy.
apply (P_approx_has_dim f).
apply has_dim_Rn.
apply gather_lm_compat, (Sigma_lm f).
apply (unisolvence_inj f).
Qed.
Lemma FE_ext :
∀ {fe1 fe2 : FE} (Hd : d fe2 = d fe1)
(Hn : ndof fe2 = ndof fe1) (Hs : shape fe2 = shape fe1),
vtx fe1 = castT (nos_eq Hd Hs) Hd (vtx fe2) →
P_approx fe1 = cast2F_fun Hd (P_approx fe2) →
Sigma fe1 = castF Hn (mapF (cast2F_fun Hd) (Sigma fe2)) →
fe1 = fe2.
Proof.
intros fe1 fe2 Hd Hn Hs; destruct fe1, fe2; simpl in *; subst.
rewrite castT_id !cast2F_fun_id castF_id; intros; subst.
f_equal; apply proof_irrel.
Qed.
Variable fe: FE.
Lemma nvtx_pos : (0 < nvtx fe)%coq_nat.
Proof.
unfold nvtx, nvtx_of_shape; case shape.
auto with arith.
apply pow_pos; auto.
Qed.
Lemma shape_dec :
{shape fe = Simplex} + {shape fe = Quad}.
Proof.
case (shape fe); [left | right]; easy.
Qed.
Lemma nvtx_dec :
{nvtx fe = S (d fe)} + {nvtx fe = (Nat.pow 2 (d fe))%nat }.
Proof.
unfold nvtx; case (shape fe); [left|right]; easy.
Qed.
Lemma nvtx_Simplex :
shape fe = Simplex → nvtx fe = S (d fe).
Proof.
unfold nvtx; case (shape fe); auto.
intros H; discriminate.
Qed.
Lemma nvtx_Quad :
shape fe = Quad → nvtx fe = (Nat.pow 2 (d fe))%nat.
Proof.
unfold nvtx; case (shape fe); auto.
intros H; discriminate.
Qed.
Definition is_local_shape_fun :=
fun (theta : '(FRd (d fe))^(ndof fe)) ⇒
(∀ i : 'I_(ndof fe), P_approx fe (theta i)) ∧
(∀ i j : 'I_(ndof fe), Sigma fe i (theta j) = kronR i j).
Lemma is_local_shape_fun_inj : ∀ theta1 theta2,
is_local_shape_fun theta1 → is_local_shape_fun theta2 →
theta1 = theta2.
Proof.
intros theta1 theta2 [H1 H1'] [H2 H2'].
pose (PC:= (has_dim_cms _ (P_approx_has_dim fe))).
generalize (lmS_bijS_val_full_equiv _ (P_approx_has_dim fe) (ndof fe) has_dim_Rn
(proj1 (gather_lm_compat _) (Sigma_lm fe)) eq_refl); intros T.
destruct T as [T1 _].
fold PC in T1.
specialize (T1 (unisolvence fe)).
apply fun_ext; intros j.
apply minus_zero_reg_sym.
generalize (cms_minus PC _ _ (H2 j) (H1 j)); intros H.
apply trans_eq with (val (mk_sub H)); try easy.
eapply trans_eq with (val zero); try easy.
apply val_eq, T1.
simpl; unfold gather, swap; simpl.
apply fun_ext; intros i.
unfold minus; rewrite (proj1 (Sigma_lm fe i)).
rewrite lm_opp.
2: apply Sigma_lm.
rewrite H1' H2'.
rewrite plus_opp_r; easy.
Qed.
d : nat ;
ndof : nat ;
d_pos : (0 < d)%coq_nat ;
ndof_pos : (0 < ndof)%coq_nat ;
shape : shape_type ;
nvtx : nat := nvtx_of_shape d shape ;
vtx : '('R^d)^nvtx ;
K_geom : 'R^d → Prop := convex_envelop vtx ;
P_approx : FRd d → Prop ;
P_approx_has_dim : has_dim P_approx ndof ;
Sigma : '(FRd d → R)^ndof ;
Sigma_lm : ∀ i, lin_map (Sigma i) ;
unisolvence_inj : KerS0 P_approx (gather Sigma) ;
}.
Lemma unisolvence : ∀ (f:FE),
bijS (P_approx f) fullset (gather (Sigma f)).
Proof.
intros f.
apply lmS_bijS_KerS0 with (ndof f) (ndof f); try easy.
apply (P_approx_has_dim f).
apply has_dim_Rn.
apply gather_lm_compat, (Sigma_lm f).
apply (unisolvence_inj f).
Qed.
Lemma FE_ext :
∀ {fe1 fe2 : FE} (Hd : d fe2 = d fe1)
(Hn : ndof fe2 = ndof fe1) (Hs : shape fe2 = shape fe1),
vtx fe1 = castT (nos_eq Hd Hs) Hd (vtx fe2) →
P_approx fe1 = cast2F_fun Hd (P_approx fe2) →
Sigma fe1 = castF Hn (mapF (cast2F_fun Hd) (Sigma fe2)) →
fe1 = fe2.
Proof.
intros fe1 fe2 Hd Hn Hs; destruct fe1, fe2; simpl in *; subst.
rewrite castT_id !cast2F_fun_id castF_id; intros; subst.
f_equal; apply proof_irrel.
Qed.
Variable fe: FE.
Lemma nvtx_pos : (0 < nvtx fe)%coq_nat.
Proof.
unfold nvtx, nvtx_of_shape; case shape.
auto with arith.
apply pow_pos; auto.
Qed.
Lemma shape_dec :
{shape fe = Simplex} + {shape fe = Quad}.
Proof.
case (shape fe); [left | right]; easy.
Qed.
Lemma nvtx_dec :
{nvtx fe = S (d fe)} + {nvtx fe = (Nat.pow 2 (d fe))%nat }.
Proof.
unfold nvtx; case (shape fe); [left|right]; easy.
Qed.
Lemma nvtx_Simplex :
shape fe = Simplex → nvtx fe = S (d fe).
Proof.
unfold nvtx; case (shape fe); auto.
intros H; discriminate.
Qed.
Lemma nvtx_Quad :
shape fe = Quad → nvtx fe = (Nat.pow 2 (d fe))%nat.
Proof.
unfold nvtx; case (shape fe); auto.
intros H; discriminate.
Qed.
Definition is_local_shape_fun :=
fun (theta : '(FRd (d fe))^(ndof fe)) ⇒
(∀ i : 'I_(ndof fe), P_approx fe (theta i)) ∧
(∀ i j : 'I_(ndof fe), Sigma fe i (theta j) = kronR i j).
Lemma is_local_shape_fun_inj : ∀ theta1 theta2,
is_local_shape_fun theta1 → is_local_shape_fun theta2 →
theta1 = theta2.
Proof.
intros theta1 theta2 [H1 H1'] [H2 H2'].
pose (PC:= (has_dim_cms _ (P_approx_has_dim fe))).
generalize (lmS_bijS_val_full_equiv _ (P_approx_has_dim fe) (ndof fe) has_dim_Rn
(proj1 (gather_lm_compat _) (Sigma_lm fe)) eq_refl); intros T.
destruct T as [T1 _].
fold PC in T1.
specialize (T1 (unisolvence fe)).
apply fun_ext; intros j.
apply minus_zero_reg_sym.
generalize (cms_minus PC _ _ (H2 j) (H1 j)); intros H.
apply trans_eq with (val (mk_sub H)); try easy.
eapply trans_eq with (val zero); try easy.
apply val_eq, T1.
simpl; unfold gather, swap; simpl.
apply fun_ext; intros i.
unfold minus; rewrite (proj1 (Sigma_lm fe i)).
rewrite lm_opp.
2: apply Sigma_lm.
rewrite H1' H2'.
rewrite plus_opp_r; easy.
Qed.
Definition shape_fun : '(FRd (d fe))^(ndof fe) := predual_basis (unisolvence fe).
Lemma shape_fun_correct : is_local_shape_fun shape_fun.
Proof.
split.
apply predual_basis_in_sub.
intros; apply predual_basis_dualF.
Qed.
Lemma shape_fun_basis : basis (P_approx fe) shape_fun.
Proof.
apply predual_basis_basis.
apply fe.
apply (Sigma_lm fe).
Qed.
Lemma shape_fun_correct : is_local_shape_fun shape_fun.
Proof.
split.
apply predual_basis_in_sub.
intros; apply predual_basis_dualF.
Qed.
Lemma shape_fun_basis : basis (P_approx fe) shape_fun.
Proof.
apply predual_basis_basis.
apply fe.
apply (Sigma_lm fe).
Qed.
[Ern]
Def 4.4, p. 78.
Definition local_interp : FRd (d fe) → FRd (d fe) :=
fun v ⇒ lin_comb (fun i ⇒ Sigma fe i v) shape_fun.
Lemma local_interp_is_poly : ∀ v : FRd (d fe),
(P_approx fe) (local_interp v).
Proof.
intros v; unfold local_interp.
rewrite (proj1 shape_fun_basis); easy.
Qed.
Lemma local_interp_lm :
lin_map local_interp.
Proof.
split.
intros v1 v2.
unfold local_interp.
rewrite -lc_plus_l.
apply lc_ext_l; intros i.
apply (Sigma_lm fe).
intros v l.
unfold local_interp.
rewrite -lc_scal_l.
apply lc_ext_l; intros i.
apply (Sigma_lm fe).
Qed.
Lemma local_interp_shape_fun : ∀ i : 'I_(ndof fe),
local_interp (shape_fun i) = shape_fun i.
Proof.
intros j.
unfold local_interp.
apply trans_eq with
(lin_comb (fun i : 'I_(ndof fe) ⇒ kronR i j) shape_fun).
apply lc_scalF_compat, extF; intros i.
f_equal; apply extF; intro.
apply shape_fun_correct.
apply lc_kron_l_in_l.
Qed.
fun v ⇒ lin_comb (fun i ⇒ Sigma fe i v) shape_fun.
Lemma local_interp_is_poly : ∀ v : FRd (d fe),
(P_approx fe) (local_interp v).
Proof.
intros v; unfold local_interp.
rewrite (proj1 shape_fun_basis); easy.
Qed.
Lemma local_interp_lm :
lin_map local_interp.
Proof.
split.
intros v1 v2.
unfold local_interp.
rewrite -lc_plus_l.
apply lc_ext_l; intros i.
apply (Sigma_lm fe).
intros v l.
unfold local_interp.
rewrite -lc_scal_l.
apply lc_ext_l; intros i.
apply (Sigma_lm fe).
Qed.
Lemma local_interp_shape_fun : ∀ i : 'I_(ndof fe),
local_interp (shape_fun i) = shape_fun i.
Proof.
intros j.
unfold local_interp.
apply trans_eq with
(lin_comb (fun i : 'I_(ndof fe) ⇒ kronR i j) shape_fun).
apply lc_scalF_compat, extF; intros i.
f_equal; apply extF; intro.
apply shape_fun_correct.
apply lc_kron_l_in_l.
Qed.
[Ern]
Eq (4.6), p. 78.
Lemma local_interp_proj : ∀ p : FRd (d fe),
(P_approx fe) p → local_interp p = p.
Proof.
intros p Hp.
rewrite (proj1 shape_fun_basis) in Hp.
inversion Hp as [ L HL ].
rewrite f_lc_compat_lm.
apply lc_ext_r; intros i.
apply local_interp_shape_fun.
apply local_interp_lm.
Qed.
(P_approx fe) p → local_interp p = p.
Proof.
intros p Hp.
rewrite (proj1 shape_fun_basis) in Hp.
inversion Hp as [ L HL ].
rewrite f_lc_compat_lm.
apply lc_ext_r; intros i.
apply local_interp_shape_fun.
apply local_interp_lm.
Qed.
[Ern]
p. 78.
Lemma local_interp_idem : ∀ v : FRd (d fe),
local_interp (local_interp v) = local_interp v.
Proof.
intros v; rewrite local_interp_proj; try easy.
apply local_interp_is_poly; try easy.
Qed.
local_interp (local_interp v) = local_interp v.
Proof.
intros v; rewrite local_interp_proj; try easy.
apply local_interp_is_poly; try easy.
Qed.
[Ern]
Eq (4.4), p. 77.
Lemma Sigma_local_interp : ∀ (i : 'I_(ndof fe))
(v : FRd (d fe)),
Sigma fe i (local_interp v) = Sigma fe i v.
Proof.
intros i v.
unfold local_interp.
rewrite f_lc_compat_lm.
apply trans_eq with
(lin_comb (fun j ⇒ kronR i j)(fun j ⇒ Sigma fe j v)).
apply lc_scalF_compat; rewrite scalF_comm_R.
f_equal; apply extF; intro; apply shape_fun_correct.
rewrite lc_kron_l_in_r; easy.
apply (Sigma_lm fe).
Qed.
End FE_Facts1.
Section FE_Facts2.
Lemma shape_fun_ext :
∀ {fe1 fe2 : FE} (H : fe1 = fe2),
shape_fun fe1 =
castF (eq_sym (f_equal ndof H))
(mapF (castF_fun (eq_sym (f_equal d H))) (shape_fun fe2)).
Proof. intros; subst; rewrite castF_id castF_fun_id; easy. Qed.
End FE_Facts2.
(v : FRd (d fe)),
Sigma fe i (local_interp v) = Sigma fe i v.
Proof.
intros i v.
unfold local_interp.
rewrite f_lc_compat_lm.
apply trans_eq with
(lin_comb (fun j ⇒ kronR i j)(fun j ⇒ Sigma fe j v)).
apply lc_scalF_compat; rewrite scalF_comm_R.
f_equal; apply extF; intro; apply shape_fun_correct.
rewrite lc_kron_l_in_r; easy.
apply (Sigma_lm fe).
Qed.
End FE_Facts1.
Section FE_Facts2.
Lemma shape_fun_ext :
∀ {fe1 fe2 : FE} (H : fe1 = fe2),
shape_fun fe1 =
castF (eq_sym (f_equal ndof H))
(mapF (castF_fun (eq_sym (f_equal d H))) (shape_fun fe2)).
Proof. intros; subst; rewrite castF_id castF_fun_id; easy. Qed.
End FE_Facts2.