Asymptotic Exactness Of The Least-Squares Finite Element .

3y ago
24 Views
2 Downloads
645.46 KB
21 Pages
Last View : 13d ago
Last Download : 3m ago
Upload by : Abby Duckworth
Transcription

SIAM J. NUMER. ANAL.Vol. 56, No. 4, pp. 2008--2028c 2018 Society for Industrial and Applied Mathematics\bigcircASYMPTOTIC EXACTNESS OF THE LEAST-SQUARESFINITE ELEMENT RESIDUAL\astCARSTEN CARSTENSEN\dagger AND JOHANNES STORN\ddaggerAbstract. The discrete minimal least-squares functional LS(f ; U ) is equivalent to the squarederror u - U 2 in least-squares finite element methods and so leads to an embedded reliable and efficient a posteriori error control. This paper enfolds a spectral analysis to prove that this natural errorestimator is asymptotically exact in the sense that the ratio LS(f ; U )/ u - U 2 tends to one as theunderlying mesh-size tends to zero for the Poisson model problem, the Helmholtz equation, the linear elasticity, and the time-harmonic Maxwell equations with all kinds of conforming discretizations.Some knowledge about the continuous and the discrete eigenspectrum allows for the computationof a guaranteed error bound C(\scrT )LS(f ; U ) with a reliability constant C(\scrT ) \leq 1/\alpha smaller thanthat from the coercivity constant \alpha . Numerical examples confirm the estimates and illustrate theperformance of the novel guaranteed error bounds with improved efficiency.Key words. least-squares finite element method, global upper bound, asymptotically exact errorestimation, sharpened reliability constants, spectral analysis, Poisson model problem, Helmholtzequation, linear elasticity, Maxwell equationsAMS subject classifications. 65N12, 65N15, 65N30DOI. 10.1137/17M11259721. Introduction. The least-squares finite element method (LSFEM) approximates the exact solution u \in X to a partial differential equation by the discreteminimizer U \in X(\scrT ) of a least-squares functional LS(f ; \bullet ) over a discrete subspaceX(\scrT ) \subset X. For the problems in this paper, namely the Poisson model problem, theHelmholtz equation, the linear elasticity, and the Maxwell equations, the functionalLS(f ; \bullet ) is equivalent to the norm \ \bullet \ 2X in X with equivalence constants \alpha and \beta . Inparticular, the discrete minimizer U \in X(\scrT ) satisfies \alpha \leq LS(f ; U )/\ u - U \ 2X \leq \betaand the computable residual LS(f ; U ) leads to a guaranteed upper bound (GUB)\ u - U \ 2X \leq \alpha - 1 LS(f ; U ) [3]. Table 1 displays computed upper and lower bounds ofthe quotient LS(f ; U )/\ u - U \ 2X for a Poisson model problem and provides numericalevidence of asymptotic exactness of the least-squares residual LS(f ; U ). This experiment suggests that the GUB \alpha - 1 LS(f ; U ) is too pessimistic for \alpha - 1 1.442114.The first main result of this paper verifies that the ratio LS(f ; U )/\ u - U \ 2Xwith the unique exact (resp., discrete) minimizer u (resp., U ) tends to one in themodel problems from section 2 as the maximal mesh size \delta of the underlying regular\ast Received by the editors April 17, 2017; accepted for publication (in revised form) April 30, 2018;published electronically July 3, 97.htmlFunding: The research of the first author was supported by the Deutsche Forschungsgemeinschaft in the Priority Program 1748 Reliable simulation techniques in solid mechanics. Developmentof non-standard discretization methods, mechanical and mathematical analysis"" under the project Foundation and application of generalized mixed FEM towards nonlinear problems in solid mechanics"" (CA 151/22-1). The research of the second author was supported by the Studienstiftung desdeutschen Volkes.\dagger Corresponding author. Institut f\"ur Mathematik, Humboldt-Universit\"at zu Berlin, D-10099Berlin, Germany (cc@math.hu-berlin.de).\ddagger Institut f\"ur Mathematik, Humboldt-Universit\"at zu Berlin, D-10099 Berlin, Germany (storn@math.hu-berlin.de).2008

2009ASYMPTOTIC EXACTNESS OF LS RESIDUALSTable 1Guaranteed lower bounds (LB) and upper bounds (UB) for the quotient LS(f ; U )/\ u - U \ 2X inthe Poisson model problem with right-hand side f \equiv 1 on the L-shaped domain \Omega ( - 1, 1)2 \setminus [0, 1)2from subsection .001112471.00070662triangulation \scrT , written \scrT \in \BbbT (\delta ), tends to zero:(1) \forall \varepsilon 0 \exists \delta 0 \forall \scrT \in \BbbT (\delta )(1 - \varepsilon ) \ u - U \ 2X \leq LS(f ; U ) \leq (1 \varepsilon ) \ u - U \ 2X .One key observation is that \varepsilon and \delta are independent of the right-hand side f in L2 (\Omega )and do not depend on the polynomial degrees of a balanced or unbalanced conforming discretization (but certainly depend on the domain and the parameters in thedifferential operators). To the best of the authors' knowledge, this is the first resultof the asymptotically exact error estimation for those problems with standard discretizations; the results in [8] are caused by an unbalanced discretization. The proofof (1) in section 3 utilizes a spectral decomposition of the ansatz space X and theGalerkin orthogonality of the error u - U . The asymptotic exactness result impliesthe overestimation of \ u - U \ 2X by the natural GUB \alpha - 1 LS(f ; U ) with the factor\alpha - 1 1 as the maximal mesh-size tends to zero. The second aim of this paper is toovercome this inefficiency by an (offline) improvement of the reliability constant C(\scrT )with \ u - U \ 2X \leq C(\scrT )LS(f ; U ) in a GUB (displayed in Figure 1), which capturesthe convergence of the least-squares residual to the exact error. Section 4 combines apriori knowledge of the continuous eigenspectrum with additional information on thediscrete eigenspectrum and achieves a computable constant C(\scrT ). The proof utilizesthe Galerkin orthogonality of the discrete solution U and so the GUB requires anexact solve but is independent of the data f ; i.e., the constant C(\scrT ) depends onlyon \scrT and C(\scrT \ ) \leq C(\scrT ) for any refinement \scrT \ of \scrT even with polynomial enrichmentof the discrete ansatz space X(\scrT \ ). A three-stage algorithm leads in subsection 4.2to C(\scrT ) and a significant improvement of the GUB \alpha - 1 LS(f ; U ), which is up to 132times larger than C(\scrT )LS(f ; U ) in Figure 1. Further numerical experiments in section 5 on the Laplace, Helmholtz, and Maxwell equations investigate the improvementin computational benchmarks: Once the relevant eigenfunctions of the least-squaressystem are resolved with sufficient accuracy, the novel reliability constant C(\scrT ) leadsto a significant improvement of the GUB. The relevant eigenmodes are of low frequency in the Poisson and elasticity problems, while certain parameters of \omega in theHelmholtz and Maxwell equations might lead to relevant high-frequency eigenmodessolely resolved for very fine meshes.Standard notation on Lebesque and Sobolev spaces applies throughout this paper,H 1 (\Omega ) : \{ v \in L2 (\Omega ; \BbbR ) : \nabla v \in L2 (\Omega ; \BbbR d )\} , H01 (\Omega ) : \{ v \in H 1 (\Omega ) : v \partial \Omega 0\} ,H(div, \Omega ) : \{ q \in L2 (\Omega ; \BbbR d ) : div q \in L2 (\Omega ; \BbbR )\} , and, for d 3 only, H(curl, \Omega ) : \{ F \in L2 (\Omega ; \BbbR 3 ) : curl F \in L2 (\Omega ; \BbbR 3 )\} , H0 (curl, \Omega ) : \{ F \in H(curl, \Omega ) : \nu \times F 0 on\partial \Omega \} with outer unit normal vector \nu \in \BbbR 3 .2. Four applications of the LSFEM. This section introduces the model problems and their finite element discretizations for a bounded polyhedral Lipschitz do-

2010CARSTEN CARSTENSEN AND JOHANNES STORN1032ku U kXLS(f ; U )α 1 LS(f ; U )C(T )LS(f ; U )10 110 5101102103104ndof105106107Fig. 1. Convergence history plot from subsection 5.4 of the squared error, residual, and GUBfor the Helmholtz equation with \omega 4.main \Omega \subset \BbbR d .2.1. Poisson model problem. Given f \in L2 (\Omega ; \BbbR ), the Poisson model problemseeks (u, p) \in X : H01 (\Omega ) \times H(div, \Omega ) with- div p f in \Omegaand\nabla u p in \Omega .First-order systems least-squares (FOSLS) methods such as, e.g., those in [2, 10, 18,20] utilize the equivalence of the Poisson model problem to the minimization of theleast-squares functionalLS(f ; v, q) : \ q - \nabla v\ 2L2 (\Omega ) \ f div q\ 2L2 (\Omega )over all (v, q) \in X with norm \ (v, q)\ 2X : \ \nabla v\ 2L2 (\Omega ) \ q\ 2L2 (\Omega ) \ div q\ 2L2 (\Omega ) .2.2. Helmholtz equation. Given some f \in L2 (\Omega ; \BbbR ) and a frequency \omega 2 0different from a Dirichlet eigenvalue of the Laplace operator, the Helmholtz equationseeks (u, p) \in X : H01 (\Omega ) \times H(div, \Omega ) with- div p - \omega 2 u f in \Omegaand\nabla u p in \Omega .This problem is well posed. The equivalent FOSLS formulation from [10] minimizesthe least-squares functionalLS(f ; v, q) : \ q - \nabla v\ 2L2 (\Omega ) \ f \omega 2 v div q\ 2L2 (\Omega )over all (v, q) \in X with norm as in subsection 2.1.2.3. Linear elasticity. Given f \in L2 (\Omega ; \BbbR d ), the linear elasticity seeks thesolution (u, \sigma ) \in X : H01 (\Omega ; \BbbR d ) \times H(div, \Omega ; \BbbR d\times d ) to- div \sigma fand\sigma \BbbC \varepsilon (u)with the linear Green strain tensor \varepsilon (u) : (\nabla u (\nabla u)\top )/2, positive Lam\'e constants\lambda and \mu , and the fourth-order elasticity tensor \BbbC [14]. The problem is equivalent tothe minimization of the least-squares functionalLS(f ; v, \tau ) : \ \BbbC - 1/2 \tau - \BbbC 1/2 \varepsilon (v)\ 2L2 (\Omega ) \ f div \tau \ 2L2 (\Omega )over all (v, \tau ) \in X with norm \ (v, \tau )\ 2X : \ \BbbC 1/2 \varepsilon (v)\ 2L2 (\Omega ) \ \BbbC - 1/2 \tau \ 2L2 (\Omega ) \ div \tau \ 2L2 (\Omega ) [9].

ASYMPTOTIC EXACTNESS OF LS RESIDUALS20112.4. Time-harmonic Maxwell equations. Given some right-hand side f \inL2 (\Omega ; \BbbR 3 ) and a frequency \omega 2 0 different from an eigenvalue of the resonant cavity problem, the time-harmonic Maxwell equations in d 3 space dimensions seek(E, H) \in X : H0 (curl, \Omega ) \times H(curl, \Omega ) with- \omega 2 E curl H f in \Omegaandcurl E - H 0 in \Omega .The problem is well posed and its solution minimizes the least-squares functionalLS(f ; F, G) : \ G - curl F \ 2L2 (\Omega ) \ f \omega 2 F - curl G\ 2L2 (\Omega )over all (F, G) \in X with norm \ (F, G)\ 2X : \omega 4 \ F \ 2L2 (\Omega ) \ curl F \ 2L2 (\Omega ) \ G\ 2L2 (\Omega ) \ curl G\ 2L2 (\Omega ) . This problem is related to the problem in [6] with the exception of anadditional term similar to the extra term in subsection 3.3.2.5. Discretization. Let \BbbT be the set of admissible and shape regular triangulations of the polyhedral bounded Lipschitz domain \Omega \subset \BbbR d into simplices[5, Chap. 5]. Given \delta 0, the subset \BbbT (\delta ) \subset \BbbT consists of all triangulations\scrT \in \BbbT with diameter hT : diam(T ) \delta for all T \in \scrT . Let \BbbP k (T ; \BbbR \ell ) denote the set of polynomials of total degree at most k \in \BbbN 0 seen as a map fromT to \BbbR \ell , \ell \in \BbbN , and define RTk (T ) : \BbbP k (T ; \BbbR \ell ) \BbbP k (T ; \BbbR ) id \subset \BbbP k (T ; \BbbR \ell ) and\scrN k (T ) : \BbbP k (T ; \BbbR 3 ) \BbbP k (T ; \BbbR 3 ) \times id \subset \BbbP k (T ; \BbbR 3 ) with the identity id on T . Definefor all k \in \BbbN 0 the Courant, Raviart--Thomas, and N\'ed\'elec element spacesS k 1 (\scrT ) : \{ V \in H 1 (\Omega ) : \forall T \in \scrT , V T \in \BbbP k 1 (T ; \BbbR )\} ,RTk (\scrT ) : \{ Q \in H(div, \Omega ) : \forall T \in \scrT , Q T \in RTk (T )\} ,\scrN k (\scrT ) : \{ F \in H(curl, \Omega ) : \forall T \in \scrT , F T \in \scrN k (T )\} .Furthermore, set S0k 1 (\scrT ) : S k 1 (\scrT ) \cap H01 (\Omega ) and \scrN 0k (\scrT ) : \scrN k (\scrT ) \cap H0 (curl, \Omega ).It is well known and understood throughout this paper that the discrete spaces X(\scrT )in Table 2 and the continuous spaces X satisfy the pointwise density property [4, 5, 7](D)\forall \varepsilon 0 \forall w \in X \exists \delta 0 \forall \scrT \in \BbbT (\delta ) \exists W \in X(\scrT )\ w - W \ X \varepsilon .3. Proof of the asymptotic exactness. The unifying analysis departs withan abstract framework and thereafter applies it to the model examples of section 2.3.1. An abstract setting. This subsection provides an abstract asymptoticexactness result based on three hypotheses.(H1) Suppose a : X \times X \rightarrow \BbbR is a scalar product that is equivalent to the scalarproduct b on the real Hilbert space (X, b) with associated norm \ \bullet \ b \ \bullet \ X . Inparticular, there exist positive constants \alpha , \beta with(2)\forall x \in X\alpha \ x\ 2b \leq a(x, x) : \ x\ 2a \leq \beta \ x\ 2b .(H2) Suppose that there exist countably many pairwise distinct positive numbers\mu (0) 1, \mu (1), \mu (2), \mu (3), . . . with closed eigenspaces E(\mu (j)) \subset X for j \in \BbbN 0 and(3)\forall j \in \BbbN 0 \forall \phi j \in E(\mu (j)) \forall x \in Xa(\phi j , x) \mu (j)b(\phi j , x).Let the eigenspaces have finite dimension dim E(\mu (j)) \in \BbbN for all j \in \BbbN (whiledim E(\mu (0)) \in \BbbN 0 \cup \{ \infty \} may be infinity or zero), and suppose that the linear hull ofall eigenspaces E(\mu (0)), E(\mu (1)), . . . is dense in X,(4)X span\{ E(\mu (j)) : j \in \BbbN 0 \} .

20121.CARSTEN CARSTENSEN AND JOHANNES STORN(H3) Suppose \mu (0) 1 is the only accumulation point of (\mu (j))j\in \BbbN 0 , limj\rightarrow \infty \mu (j) Given a right-hand side F \in X \ast in the dual X \ast of X, let u \in X be the uniquesolution to a(u, v) F (v) for all v \in X. Furthermore, let X(\scrT ) \subset X satisfy thedensity property (D), and define the discrete solution U \in X(\scrT ) with a(U, V ) F (V )for all V \in X(\scrT ).Theorem 3.1. Suppose (H1)--(H3), (D), and \varepsilon 0. Then there exists some\delta 0 for all F \in X \ast such that for all \scrT \in \BbbT (\delta )(5)(1 - \varepsilon )\ u - U \ 2b \leq \ u - U \ 2a \leq (1 \varepsilon )\ u - U \ 2b .Some remarks are in order before the proof of the theorem concludes this subsection.Remark 3.2 (bounded eigenvalues). It follows from (H1) that \mu (0), \mu (1), \mu (2), . . .are bounded in the compact interval [\alpha , \beta ] and, since \mu (0) 1, it holds that \alpha \leq 1 \leq \beta .Remark 3.3 (orthogonal eigenspaces). The eigenvectors \phi j \in E(\mu (j)) and \phi k \inE(\mu (k)) with j, k \in \BbbN satisfy \mu (j)b(\phi j , \phi k ) a(\phi j , \phi k ) a(\phi k , \phi j ) \mu (k)b(\phi k , \phi j ).If j \not k, it holds that 0 \not \mu (j) \not \mu (k) \not 0, and so b(\phi j , \phi k ) a(\phi j , \phi k ) 0. Thus,(6)\forall j, k \in \BbbN 0 \wedge j \not kE(\mu (j)) \bot a E(\mu (k)) and E(\mu (j)) \bot b E(\mu (k)).Remark 3.4 (orthogonal decomposition of X). Given an index set J \subset \BbbN 0 , defineX(J) as the closure of span\{ E(\mu (j)) : j \in J\} , and set the complement J c : \BbbN 0 \setminus J.Then\sum(4) implies that any v \in X can\sum be decomposed into v w z with somew j\in J wj \in X(J) and some z k\in J c zk \in X(J c ) such that wj \in E(\mu (j)) forall j \in J and zk \in E(\mu (k)) for all k \in J c . Since X(J) and X(J c ) are closed withrespect to the norm \ \bullet \ b , (6) implies b(w, z) 0. This proves the b-orthogonalityX(J) \bot b X(J c ). Similar arguments and the equivalence (2) of \ \bullet \ b and \ \bullet \ a implythe a-orthogonality X(J) \bot a X(J c ).Remark 3.5 (built-in error control of LSFEMs). The least-squares formulationsfrom section 2 allow (H1)--(H3) such that \ u - U \ 2a LS(f ; U ) is a computableresidual and serves as an error estimator for the unknown error \ u - U \ b \ u - U \ X .The ellipticity in (H1) leads to(7)\alpha \ u - U \ 2b \leq LS(f ; U ) \ u - U \ 2a \leq \beta \ u - U \ 2b .This is well known in the least-squares community and called reliability and efficiencyin the a posteriori error analysis. It is a consequence of Theorem 3.1 that the GUBin (7) leads to an overestimation by the factor \alpha - 1 as the mesh size tends to zero.Proof of Theorem 3.1. The Galerkin orthogonality a(u - U, W ) 0 for all W \inX(\scrT ) is rewritten as u - U \in X(\scrT )\bot : \{ v \in X : \forall W \in X(\scrT ), a(v, W ) 0\} . Thenthe theorem follows from the more general assertion(8)\forall \varepsilon 0 \exists \delta 0 \forall \scrT \in \BbbT (\delta ) \forall v \in X(\scrT )\bot(1 - \varepsilon )\ v\ 2b \leq \ v\ 2a \leq (1 \varepsilon )\ v\ 2b .To prove (8), let 0 \varepsilon 1 and v \in X(\scrT )\bot with \ v\ b 1.Step 1 (decomposition of v). Recall \mu (0), \mu (1), . . . from (H2), and, given \varepsilon 0,define the index set J(\varepsilon ) : \{ j \in \BbbN : 1 - \mu (j) \varepsilon \} with complement J c (\varepsilon ) : \BbbN 0 \setminus J(\varepsilon ). It is a consequence of (H3) that the index set J(\varepsilon ) is finite. As outlined

2013ASYMPTOTIC EXACTNESS OF LS RESIDUALSTable 2Notation in subsection 3.2.PoissonHelmholtzElasticityMaxwell\BbbM\BbbR d\BbbR d\BbbR d\times d\BbbR 3Aidid\BbbCid\gamma0\omega 20\omega 2D\ast- div- div- divcurlD\nabla\nabla\varepsilon (\bullet )curlX(\scrT )S0k 1 (\scrT ) \times RTk (\scrT )S0k 1 (\scrT ) \times RTk (\scrT )S0k 1 (\scrT )d \times RTk (\scrT )d\scrN 0k (\scrT ) \times \scrN k (\scrT )in Remark 3.4, (H2) leads to the a- and b-orthogonal decomposition v w z withw \in X(J(\varepsilon )) and z \in X(J c (\varepsilon )). The Pythagoras theorem reads1 \ v\ 2b \ w\ 2b \ z\ 2b(9)\ v\ 2a \ w\ 2a \ z\ 2a .andStep 2 (upper bound for \ w\ a ). Let (\phi 1 , . . . , \phi m ) \sumbe a b-orthonormal basis ofmspan\{ E(\mu (j)) : j \in J(\varepsilon )\} span\{ \phi 1 , . . . , \phi m \} with w k 1 \xi k \phi k . The density (D)leads to \delta 0 such that\scrT \in \BbbT (\delta ) there exists a \Phi k \in X(\scrT )\sum m\surd for all k 1, . . . , m andwith \ \phi k - \Phi k \ b \leq \varepsilon / m. The discrete W : k 1 \xi k \Phi k \in X(\scrT ) satisfies\ w - W \ b \leqm\sumk 1- 1/2 \xi k \ \phi k - \Phi k \ b \leq m\varepsilonm\sum \xi k \leq \varepsilon\biggl( \summk 1\xi k2\biggr) 1/2 \varepsilon \ w\ b .k 1The combination with a(w, z) 0 a(v, W ), a Cauchy--Schwarz inequality, and (2)proves\ w\ 2a a(w, v) a(w - W, v) \leq \beta \ w - W \ b \leq \varepsilon \beta \ w\ b \leq \alpha - 1/2 \varepsilon \beta \ w\ a .(10)Step 3 (upper and lower bounds for \ z\ 2a ). Since z is in the closure of the linearchull: j \in J\sum(\varepsilon )\} with respect to \ \bullet \ a and \ \bullet \ b , the sums \ z\ 2a \sum span\{ E(\mu (j))222j\in J c (\varepsilon ) \ zj \ a and \ z\ b j\in J c (\varepsilon ) \ zj \ b converge. Then 1 - \varepsilon \leq \mu (j) \leq 1 \varepsilon and\ zj \ 2a \mu (j)\ zj \ 2b for all j \in J c (\varepsilon ) imply(1 - \varepsilon )\ z\ 2b \leq \ z\ 2a \leq (1 \varepsilon )\ z\ 2b .(11)Step 4 (upper bound for \ v\ 2a ). The combination of (9)--(11) proves\ v\ 2a \ z\ 2a \ w\ 2a \leq (1 \varepsilon )\ z\ 2b \ w\ 2a \leq 1 \varepsilon \varepsilon 2 \beta 2 /\alpha .St

The least-squares finite element method (LSFEM) approxi-mates the exact solution u \in X to a partial differential equation by the discrete minimizer U\in X(\scrT ) of a least-squares functional LS(f;\bullet ) over a discrete subspace X(\scrT ) \subset X. For the problems in this paper, namely the Poisson model problem, the

Related Documents:

May 02, 2018 · D. Program Evaluation ͟The organization has provided a description of the framework for how each program will be evaluated. The framework should include all the elements below: ͟The evaluation methods are cost-effective for the organization ͟Quantitative and qualitative data is being collected (at Basics tier, data collection must have begun)

Silat is a combative art of self-defense and survival rooted from Matay archipelago. It was traced at thé early of Langkasuka Kingdom (2nd century CE) till thé reign of Melaka (Malaysia) Sultanate era (13th century). Silat has now evolved to become part of social culture and tradition with thé appearance of a fine physical and spiritual .

On an exceptional basis, Member States may request UNESCO to provide thé candidates with access to thé platform so they can complète thé form by themselves. Thèse requests must be addressed to esd rize unesco. or by 15 A ril 2021 UNESCO will provide thé nomineewith accessto thé platform via their émail address.

̶The leading indicator of employee engagement is based on the quality of the relationship between employee and supervisor Empower your managers! ̶Help them understand the impact on the organization ̶Share important changes, plan options, tasks, and deadlines ̶Provide key messages and talking points ̶Prepare them to answer employee questions

Dr. Sunita Bharatwal** Dr. Pawan Garga*** Abstract Customer satisfaction is derived from thè functionalities and values, a product or Service can provide. The current study aims to segregate thè dimensions of ordine Service quality and gather insights on its impact on web shopping. The trends of purchases have

Chính Văn.- Còn đức Thế tôn thì tuệ giác cực kỳ trong sạch 8: hiện hành bất nhị 9, đạt đến vô tướng 10, đứng vào chỗ đứng của các đức Thế tôn 11, thể hiện tính bình đẳng của các Ngài, đến chỗ không còn chướng ngại 12, giáo pháp không thể khuynh đảo, tâm thức không bị cản trở, cái được

Le genou de Lucy. Odile Jacob. 1999. Coppens Y. Pré-textes. L’homme préhistorique en morceaux. Eds Odile Jacob. 2011. Costentin J., Delaveau P. Café, thé, chocolat, les bons effets sur le cerveau et pour le corps. Editions Odile Jacob. 2010. Crawford M., Marsh D. The driving force : food in human evolution and the future.

Le genou de Lucy. Odile Jacob. 1999. Coppens Y. Pré-textes. L’homme préhistorique en morceaux. Eds Odile Jacob. 2011. Costentin J., Delaveau P. Café, thé, chocolat, les bons effets sur le cerveau et pour le corps. Editions Odile Jacob. 2010. 3 Crawford M., Marsh D. The driving force : food in human evolution and the future.