A Matlab implementation of two dimensional homogenization formulas

Yves Capdeboscq - Michael Vogelius

Introduction

This Document and the M-files mentioned within can be found

here in tgz format (or in zip format)

All files described in part A can be run with Matlab 5 or 6, without the PDE Toolbox. To run the demos listed in part B, the PDE Toolbox is required.


The implementation examples available in the Examples directory require the Matlab PDE Toolbox. The results of these computations are in the directory Examples\( \backslash \)Results.


These informations are also available as line commands: type help <function> to obtain help.



A Haar and Multiwavelet homogenization files


A.1 AnalyticWaveHom.m

Usage : [RH,RMW,RT]=WaveHom(M,N)
   
Input : \( M \) and \( N \), \( 2 \) by \( 2 \) matrices
   
  \( M \) the \( (x,x) \) periodic diffusion coefficient
   
  \( N \) the \( (y,y) \) periodic diffusion coefficient
   
  \( \frac{d}{dx}M\frac{d}{dx}+\frac{d}{dy}N\frac{d}{dy} \)


For a general symmetric periodic matrix use HaarSmoothing instead.


The output of this function is the same as the output of the HaarSmoothing function, in periodic 2 by 2 case. Here the formulas are implemented in their simpler form (in terms of the \( X,Y,Z \) variables) without explicitly using wavelet coefficients.


The output is composed of 3 row column vectors


RH, a 3 row column containing the \( ((x,x),(x,y),(y,y)) \)
  homogenized coefficients according to the Haar Wavelets
   
RMW, a 3 row column containing the \( ((x,x),(x,y),(y,y)) \)
  homogenized coefficients according to the Legendre-2 Multiwavelet
   
RT, a 3 row column containing the \( ((x,x),(x,y),(y,y)) \)
  exact homogenized formulas in the case \( M=N \).


A.2 AsymptoticRationalReduction.m

Usage : RES=AsymptoticRationalReduction(O)
   
Input : a \( L\times L\times 3\times 6\times 16\times 16 \) array,
   
  constructed by MRADecompositionMatrix (Paragraph A.8).


Part of MWSmoothing (Paragraph A.9).


The output of this function is an \( L\times L\times 3 \) array corresponding to the second order coefficients of the reduced operator taking the limit of the differentiation parameter \( \theta =1/2+\epsilon \) when epsilon tends to zero.


RES(:,:,1) corresponds to the (x,x) diffusion coefficient
   
RES(:,:,2) corresponds to the (y,y) diffusion coefficient
   
RES(:,:,3) corresponds to the (x,y) diffusion coefficient


\begin{displaymath}
\frac{d}{dx}\textrm{RES}(x,y,1)\frac{d}{dx}+\frac{d}{dy}\textrm{RES}(x,y,2)\frac{d}{dy}\end{displaymath}


\begin{displaymath}
+\frac{d}{dy}\textrm{RES}(x,y3)\frac{d}{dx}+\frac{d}{dx}\textrm{RES}(x,y,3)\frac{d}{dy}\end{displaymath}


See also MWSmoothing (Paragraph A.9) and MRADecompositionMatrix (Paragraph A.8).


A.3 Extension.m

Usage : Extension(sig)
   
Input : sig a \( n\times n\times 3 \) array corresponding to
   
  the \( (x,x) \), \( (y,y) \) and \( (x,y) \) diffusion coefficient
   
  \( \frac{d}{dx}\textrm{sig}(x,y,1)\frac{d}{dx}+\frac{d}{dy}\textrm{sig}(x,y,2)\frac{d}{dy} \)
   
  \( +\frac{d}{dx}\textrm{sig}(x,y,3)\frac{d}{dy}+\frac{d}{dy}\textrm{sig}(x,y,3)\frac{d}{dx} \)


The output is a \( (n+2)\times (n+2)\times 3 \) array, which is an extension of sig according to the following rule



\begin{displaymath}
\begin{array}{ccccccccccc}
& & & & & f & e & f & g & h & g\...
... & & & & & & & \\
& & & & & f & e & f & g & h & g
\end{array}\end{displaymath}


This corresponds to a locally periodic extension.


See also HaarSmoothing (Paragraph A.7), MWSmoothing (Paragraph A.9) and PeriodicExtension (Paragraph A.10).


A.4 HaarFormulaX.m

Usage : HaarFormulaX(a1x,a2x,a3x,a4x,a1y,a2y,a3y,a4y,b1,b2,b3,b4)
   
Input : Haar wavelet coefficients corresponding to
   
  \( a_{1}^{x}\ldots a_{4}^{x} \) the \( (x,x) \) diffusion coefficient
   
  \( a_{1}^{y}\ldots a_{4}^{y} \) the \( (y,y) \) diffusion coefficient
   
  \( b_{1}\ldots b_{4} \) the \( (x,y) \) diffusion coefficient


The output is the \( (x,x) \) Haar homogenized diffusion coefficient.


See also HaarSmoothing (Paragraph A.7), HaarFormulaY (Paragraph A.5) and HaarFormulaZ (Paragraph A.6). The formula is

\begin{displaymath}
((a_{1}^{x})^{3}a_{1}^{y}-(a_{1}^{y}a_{3}^{x})^{2}+(b_{2}(a_...
...2}-2(a_{3}^{x}a_{3}^{y}a_{4}^{x}+a_{2}^{x}a_{4}^{x}b_{2}\ldots \end{displaymath}


\begin{displaymath}
-2a_{3}^{x}b_{1}b_{2}+a_{4}^{x}b_{2}^{2}+a_{3}^{x}a_{4}^{x}b...
...a_{3}^{x}(a_{3}^{y}+b_{3}))b_{4}^{2}+2a_{4}^{x}b_{4}^{3}\ldots \end{displaymath}


\begin{displaymath}
+b_{4}^{4}+2a_{1}^{y}a_{3}^{x}(-(a_{3}^{x}b_{1})+a_{4}^{x}(a...
...a_{1}^{y})^{2}-(a_{3}^{y})^{2}+2a_{1}^{y}b_{1}-b_{2}^{2}\ldots \end{displaymath}


\begin{displaymath}
-2a_{3}^{y}b_{3}-b_{3}^{2}-b_{4}^{2})+a_{1}^{x}(-a_{1}^{y}((...
...x})^{2}+2a_{2}^{x}b_{2}+2b_{2}^{2})+2b_{2}(-(b_{1}b_{2})\ldots \end{displaymath}


\begin{displaymath}
+a_{4}^{x}(a_{3}^{y}+b_{3}))+2(-(a_{1}^{y}a_{4}^{x})+a_{3}^{...
...}+2b_{2})(a_{3}^{y}+b_{3}))b_{4}-2(a_{1}^{y}+b_{1})b_{4}^{2}))/\end{displaymath}


\begin{displaymath}
((a_{1}^{x})^{2}a_{1}^{y}-a_{1}^{y}(a_{2}^{x}+b_{2})^{2}+2(a...
..._{1})b_{4}^{2}+a_{1}^{x}((a_{1}^{y})^{2}+2a_{1}^{y}b_{1}\ldots \end{displaymath}


\begin{displaymath}
-(a_{3}^{y}+b_{3})^{2}-b_{4}^{2}))\end{displaymath}


A.5 HaarFormulaY.m

Usage : HaarFormulaY(a1x,a2x,a3x,a4x,a1y,a2y,a3y,a4y,b1,b2,b3,b4)
   
Input : Haar wavelet coefficients corresponding to
   
  \( a_{1}^{x}\ldots a_{4}^{x} \) the \( (x,x) \) diffusion coefficient
   
  \( a_{1}^{y}\ldots a_{4}^{y} \) the \( (y,y) \) diffusion coefficient
   
  \( b_{1}\ldots b_{4} \) the \( (x,y) \) diffusion coefficient


The output is the \( (y,y) \) Haar homogenized diffusion coefficient.


See also HaarSmoothing (Paragraph A.7), HaarFormulaX (Paragraph A.4) and HaarFormulaZ (Paragraph A.6). The formula is


\begin{displaymath}
((a_{1}^{x})^{2}(a_{1}^{y}-a_{2}^{y})(a_{1}^{y}+a_{2}^{y})+(a_{2}^{y}(a_{2}^{x}+b_{2})-b_{3}(a_{3}^{y}+b_{3}))^{2}\ldots \end{displaymath}


\begin{displaymath}
-2(a_{2}^{y}a_{4}^{y}(a_{2}^{x}+b_{2})+(a_{3}^{y}a_{4}^{y}-2...
...y}b_{1})b_{3}+a_{4}^{y}b_{3}^{2})b_{4}-(-(a_{4}^{y})^{2}\ldots \end{displaymath}


\begin{displaymath}
+2a_{2}^{y}(a_{2}^{x}+b_{2})+2b_{3}(a_{3}^{y}+b_{3}))b_{4}^{...
..._{1}^{y})^{2}((a_{2}^{x}+b_{2})^{2}+b_{3}^{2}+b_{4}^{2})\ldots \end{displaymath}


\begin{displaymath}
+2a_{1}^{y}(b_{3}(a_{4}^{y}(a_{2}^{x}+b_{2})-b_{1}b_{3})+(a_...
...(2a_{2}^{x}+a_{2}^{y}+2b_{2})b_{3})b_{4}-b_{1}b_{4}^{2})\ldots \end{displaymath}


\begin{displaymath}
+a_{1}^{x}((a_{1}^{y})^{3}+2(a_{1}^{y})^{2}b_{1}+2a_{2}^{y}(...
...{1})+a_{4}^{y}(a_{3}^{y}+b_{3})+(a_{3}^{y}+2b_{3})b_{4})\ldots \end{displaymath}


\begin{displaymath}
-a_{1}^{y}((a_{2}^{y})^{2}+(a_{3}^{y})^{2}+(a_{4}^{y})^{2}+2a_{3}^{y}b_{3}+2b_{3}^{2}+2a_{4}^{y}b_{4}+2b_{4}^{2})))/\end{displaymath}


\begin{displaymath}
((a_{1}^{x})^{2}a_{1}^{y}-a_{1}^{y}(a_{2}^{x}+b_{2})^{2}+2(a...
...4}-(a_{1}^{y}+2b_{1})b_{4}^{2}+a_{1}^{x}((a_{1}^{y})^{2}\ldots \end{displaymath}


\begin{displaymath}
+2a_{1}^{y}b_{1}-(a_{3}^{y}+b_{3})^{2}-b_{4}^{2}))\end{displaymath}


A.6 HaarFormulaZ.m

Usage : HaarFormulaZ(a1x,a2x,a3x,a4x,a1y,a2y,a3y,a4y,b1,b2,b3,b4)
   
Input : Haar wavelet coefficients corresponding to
   
  \( a_{1}^{x}\ldots a_{4}^{x} \) the \( (x,x) \) diffusion coefficient
   
  \( a_{1}^{y}\ldots a_{4}^{y} \) the \( (y,y) \) diffusion coefficient
   
  \( b_{1}\ldots b_{4} \) the \( (x,y) \) diffusion coefficient


The output is the \( (x,y) \) Haar homogenized diffusion coefficient.


See also HaarSmoothing (Paragraph A.7), HaarFormulaX (Paragraph A.4) and HaarFormulaY (Paragraph A.5). The formula is

\begin{displaymath}
(-(a_{2}^{y}a_{3}^{x}a_{3}^{y}b_{2})+a_{1}^{y}a_{3}^{x}a_{4}...
...y}b_{1}-a_{2}^{y}b_{2})+(a_{2}^{x})^{2}(-(a_{1}^{y}b_{1})\dots \end{displaymath}


\begin{displaymath}
+a_{2}^{y}b_{2})-(a_{1}^{y})^{2}a_{3}^{x}b_{3}+a_{3}^{x}(a_{...
...{3}-a_{3}^{y}b_{2}^{2}b_{3}+2a_{3}^{x}a_{3}^{y}b_{3}^{2}\ldots \end{displaymath}


\begin{displaymath}
-b_{2}^{2}b_{3}^{2}+a_{3}^{x}b_{3}^{3}+(a_{1}^{y}a_{2}^{y}a_...
...y}a_{4}^{x}b_{2}+2a_{3}^{y}b_{1}b_{2}-a_{4}^{y}b_{2}^{2}\ldots \end{displaymath}


\begin{displaymath}
-(a_{3}^{y}a_{4}^{x}+a_{3}^{x}a_{4}^{y}-2(a_{1}^{y}+2b_{1})b...
...{4}^{y})+b_{1}(a_{1}^{y}+2b_{1})+b_{2}(a_{2}^{y}+b_{2})+\ldots \end{displaymath}


\begin{displaymath}
(a_{3}^{x}+b_{3})(a_{3}^{y}+b_{3}))b_{4}^{2}+(a_{4}^{x}+a_{4...
...}^{y})^{2}b_{1}-(a_{3}^{y})^{2}b_{1}+2a_{1}^{y}b_{1}^{2}\ldots \end{displaymath}


\begin{displaymath}
-a_{1}^{y}a_{2}^{y}b_{2}+a_{3}^{y}a_{4}^{y}b_{2}-2a_{2}^{y}b...
...-2a_{3}^{y}b_{1}b_{3}+a_{4}^{y}b_{2}b_{3}-b_{1}b_{3}^{2}\ldots \end{displaymath}


\begin{displaymath}
+(a_{2}^{y}a_{3}^{x}+a_{2}^{y}a_{3}^{y}-a_{1}^{y}(a_{4}^{x}+...
...}+b_{1})b_{4}^{2})+a_{2}^{x}(a_{1}^{y}a_{3}^{x}a_{4}^{y}\ldots \end{displaymath}


\begin{displaymath}
-2a_{1}^{y}b_{1}b_{2}+a_{1}^{y}a_{4}^{x}b_{3}-a_{3}^{y}b_{2}...
...4}^{y}b_{2}+2b_{1}b_{3}+a_{1}^{y}(a_{3}^{x}+b_{3}))b_{4}\ldots \end{displaymath}


\begin{displaymath}
-b_{2}b_{4}^{2}-a_{2}^{y}(-2b_{2}^{2}+a_{3}^{x}(a_{3}^{y}+b_{3})+b_{4}(a_{4}^{x}+b_{4}))))/\end{displaymath}


\begin{displaymath}
((a_{1}^{x})^{2}a_{1}^{y}-a_{1}^{y}(a_{2}^{x}+b_{2})^{2}+2(a...
...{3})b_{4}-(a_{1}^{y}+2b_{1})b_{4}^{2}+a_{1}^{x}((a_{1}^{y})^{2}\end{displaymath}


\begin{displaymath}
+2a_{1}^{y}b_{1}-(a_{3}^{y}+b_{3})^{2}-b_{4}^{2}))\end{displaymath}


A.7 HaarSmoothing.m

Usage : HaarSmoothing(sig,EType)
   
Input : sig, a \( n\times n\times 3 \) dimensional array
   
  EType, an optional parameter.
   
  The first \( n\times n \) matrix represents the \( (x,x) \) diffusion coefficient
   
  The second \( n\times n \) matrix represents the \( (y,y) \) diffusion coefficient
   
  The third \( n\times n \) matrix represents the \( (x,y)=(y,x) \) diffusion coefficient
   
  If Etype is set to \( 2 \), PeriodicExtension (Paragraph A.10) is used
   
  instead of Extension (Paragraph A.3).


Part of Smoothing (Paragraph A.11).


The output is an \( n\times n\times 3 \) array, after an extension to an \( (n+2)\times (n+2)\times 3 \) matrix, and 2 applications of Haar homogenization formulas.


See also Extension (Paragraph A.3), PeriodicExtension (Paragraph A.10), WaveletSquare (Paragraph A.12), HaarFormulaX (Paragraph A.4), HaarFormulaY (Paragraph A.5) and HaarFormulaZ (Paragraph A.6).


A.8 MRADecompositionMatrix.m

Usage : MRADecompositionMatrix(sig)
   
Input : an \( L\times L\times 4\times 3\times \) array, where
   
  \( (i,j,k,l) \) is the kth Haar wavelet coefficient of the
   
  \( l \)-th diffusion coefficient:
   
  \( l=1 \) corresponds to (x,x)
   
  \( l=2 \) corresponds to (y,y)
   
  \( l=3 \) corresponds to (x,y)
   
  at the point \( (i,j) \).


The output is an \( L\times L\times 3\times 6\times 16\times 16 \) array, corresponding to the decomposition of the operator according to the multi-resolution analysis induced by the multiwavelet used.


\begin{displaymath}
\frac{d}{dx}a^{x}(x,y)\frac{d}{dx}+\frac{d}{dy}a^{y}(x,y)\fr...
...}(i,j) & \\
\ldots & & \ldots
\end{array}\right]
\end{array}\end{displaymath}


where \( \textrm{ER}(i,j)=\textrm{E}_{1}(i,j)+\textrm{E}_{2}(i,j)\frac{d}{dx}+\textrm{...
...\textrm{E}_{5}(i,j)\frac{d^{2}}{dy^{2}}+\textrm{E}_{6}(i,j)\frac{d^{2}}{dxdy} \).


There is a degree of freedom in the choice of the difference operator (see [1]), which will be referred to as the differentiation parameter, \( \theta \). A differentiation parameter \( \theta =1 \) corresponds to a standard forward-backward scheme, comparable to that used for the Haar basis -in that case, there is also a differentiation parameter but the result does not depend on the choice of this parameter. A differentiation parameter \( \theta =\frac{1}{2} \) corresponds to a centered difference scheme. The Gaussian reduction we perform cannot be used in that case. This procedure provided the array ER in the case \( \theta =\frac{1}{2}+\epsilon \), with the coefficients in \( \textrm{E}_{1}\ldots \textrm{E}_{6} \) limited to the second order in epsilon:

\begin{displaymath}
\textrm{E}=\textrm{E}^{0}+\epsilon \textrm{E}^{1}+\epsilon ^{2}\textrm{E}^{2}.\end{displaymath}


We therefore have ER of dimension

\begin{displaymath}
\begin{array}[b]{ccccccccccc}
\underbrace{L} & \times & \und...
...,\frac{d}{dx},\ldots \frac{d^{2}}{dxdy} & & i & & j
\end{array}\end{displaymath}


See also MWSmoothing (Paragraph A.9) and AsymptoticRationalReduction (Paragraph A.2).

  1. B. Alpert, G. Beylkin, D. Gines, L. Vozovoi, Adaptive Solution of Partial Differential Equations in Multiwavelet Bases, May 1999.


A.9 MWSmoothing.m

Usage : MWSmoothing(sig,EType)
   
Input : sig, a \( n\times n\times 3 \) dimensional array
   
  EType, an optional parameter.
   
  The first \( n\times n \) matrix represents the \( (x,x) \) diffusion coefficient
   
  The second \( n\times n \) matrix represents the \( (y,y) \) diffusion coefficient
   
  The third \( n\times n \) matrix represents the \( (x,y)=(y,x) \) diffusion coefficient
   
  If Etype is set to \( 2 \), PeriodicExtension (Paragraph A.10) is used
   
  instead of Extension (Paragraph A.3).


Part of Smoothing (Paragraph A.11).


The output is an \( n\times n\times 3 \) array, after an extension to an \( (n+2)\times (n+2)\times 3 \) matrix, and 2 numerical multiwavelet homogenization iterations.


See also Extension (Paragraph A.3), PeriodicExtension (Paragraph A.10), WaveletSquare (Paragraph A.12), MRADecompositionMatrix (Paragraph A.8) and AsymptoticRationalReduction (Paragraph A.2).


A.10 PeriodicExtension.m

Usage : PeriodicExtension(sig)
   
Input : sig a \( n\times n\times 3 \) array corresponding to
   
  the \( (x,x) \), \( (y,y) \) and \( (x,y) \) diffusion coefficients.


The output is a \( (n+2)\times (n+2)\times 3 \) array, which is an extension of sig according to the following rule

\begin{displaymath}
\begin{array}{ccccccccccc}
& & & & & L & I & J & K & L & I\...
... & & & & & & & \\
& & & & & D & A & B & C & D & A
\end{array}\end{displaymath}


See also HaarSmoothing (Paragraph A.7), MWSmoothing (Paragraph A.9) and Extension (Paragraph A.3).


A.11 Smoothing.m

Usage : resu=Smoothing(sig,Type,NbIter,BType)
   
  e.g. resu=Smoothing(A,'Haar','auto')
   
Input : sig a \( n\times n\times \left( 1,2\textrm{ or }3\right) \)array representing
   
  an isotropic diffusion operator coefficient
   
  \( \frac{d}{dx}\textrm{sig}(x,y)\frac{d}{dx}+\frac{d}{dy}\textrm{sig}(x,y)\frac{d}{dy} \)
   
  a diagonal diffusion operator coefficient
   
  \( \frac{d}{dx}\textrm{sig}(x,y,1)\frac{d}{dx}+\frac{d}{dy}\textrm{sig}(x,y,2)\frac{d}{dy} \)
   
  a symmetric diffusion operator coefficient
   
  \( \frac{d}{dx}\textrm{sig}(x,y,1)\frac{d}{dx}+\frac{d}{dy}\textrm{sig}(x,y,2)\frac{d}{dy}+ \)
   
  \( \frac{d}{dx}\textrm{sig}(x,y,3)\frac{d}{dy}+\frac{d}{dy}\textrm{sig}(x,y,3)\frac{d}{dx} \)
   
  type can be 'Haar' or 'MW'
   
  NbIter can be any number or 'auto' which means the smoothing method
   
  determines itself when to stop (variations between two steps less than \( 5\% \))
   
  BType (optional), if set to 2, sets globally periodic boundary
   
  condition instead of local ones.


The output is a smoothed image, homogenized locally according to the Haar or Multiwavelet homogenization formulas.


See also HaarSmoothing (Paragraph A.7), MWSmoothing (Paragraph A.9) and AnalyticWaveHom (Paragraph A.1).


A.12 WaveletSquare.m

Usage : WaveletSquare(sig)
   
Input : sig an \( n\times n\times 3 \) array corresponding to
   
  the \( (x,x) \), \( (y,y) \) and \( (x,y) \) diffusion coefficients.


The output is an \( (n-1)\times (n-1)\times 4\times 3 \) array corresponding to the value of the Haar Wavelet coefficients.

B Examples Files

The test examples provided are


These test consist of comparison between finite element solutions computed on the non averaged data and the finite element solutions (or the exact analytic solution when it is known) using the exact homogenized formulas and the averaging formulas described earlier. Use the help command for additional information.

About this document ...

A Matlab implementation of two dimensional homogenization formulas

This document was generated using the LaTeX2HTML translator Version 2K.1beta (1.47)

Copyright © 1993, 1994, 1995, 1996, Nikos Drakos, Computer Based Learning Unit, University of Leeds.
Copyright © 1997, 1998, 1999, Ross Moore, Mathematics Department, Macquarie University, Sydney.

The command line arguments were:
latex2html -no_subdir -split 0 -show_section_numbers /tmp/lyx_tmpdir11037pnIcG7/lyx_tmpbuf11037VBlN3N/Documentation.tex

The translation was initiated by Yves Capdeboscq on 2002-01-15


next_inactive up previous
Yves Capdeboscq 2002-01-15