hmm_stima_tri

PURPOSE ^

HMM_STIMA_TRI Computes element stiffness matrix for macro triangles.

SYNOPSIS ^

function [A_macro]=hmm_stima_tri(NMicro, epsilon, MicroElements,MacroVertices, Constraints, delta, bctype)

DESCRIPTION ^

HMM_STIMA_TRI   Computes element stiffness matrix for macro triangles.
   A_MACRO = HMM_STIMA_TRI(VERTICES, MACROQUADNODE, EPSILON) 
   computes the entry to the macro element stiffness matrix for macro
   triangles using parallelograms on the micro sampling domain and
   specified boundary conditions on the micro problem.

   NMICRO  is the number of degrees of freedom per space dimension in the
   micro domain

   EPSILON is the periodicity of the tensor a.

   MICROELEMENTS has dimension (NMicro-1)^2 x 4 and contain the node
   numbers of the triangulation of the micro sampling domain.

   MACROVERTICES contains the node coordinates of the macro triangles.

   CONSTRAINTS contains information about the micro boundary problem. The
   format depends on whether dirichlet or periodic boundary conditions are
   used. See make_contraints.m for details.

   DELTA is the width or height of the square sampling domain.

   BCTYPE is either "dirichlet" or "periodic" and determines the boundary
   condition of the micro problem.
   
   A_MACRO has dimension 4 x 4. The vertices are numbered 
   anti-clockwise. 

   This function should not be modified.


   The code is available at http://anmc.epfl.ch/ and described in 
   further detail in 

   A. Abdulle and A. Nonnenmacher
   "A short and versatile finite element multiscale code for
   homogenization problems"
   Computer Methods in Applied Mechanics and Engineering,
   http://dx.doi.org/10.1016/j.cma.2009.03.019

   Please cite this article in any publication describing research
   performed using the software.


   Email           : assyr.abdulle@epfl.ch and achim.nonnenmacher@epfl.ch
   Last updated    : 04/29/2009 with MATLAB 7.4

   FE_HMM2D is Copyright (C) 2009 A. Abdulle and A. Nonnenmacher. 
   The software is provided free for non-commercial use unter the terms of 
   the GNU General Public License. See "copyright.m" for full details.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [A_macro]=hmm_stima_tri(NMicro, epsilon, MicroElements,...
0002     MacroVertices, Constraints, delta, bctype)
0003 %HMM_STIMA_TRI   Computes element stiffness matrix for macro triangles.
0004 %   A_MACRO = HMM_STIMA_TRI(VERTICES, MACROQUADNODE, EPSILON)
0005 %   computes the entry to the macro element stiffness matrix for macro
0006 %   triangles using parallelograms on the micro sampling domain and
0007 %   specified boundary conditions on the micro problem.
0008 %
0009 %   NMICRO  is the number of degrees of freedom per space dimension in the
0010 %   micro domain
0011 %
0012 %   EPSILON is the periodicity of the tensor a.
0013 %
0014 %   MICROELEMENTS has dimension (NMicro-1)^2 x 4 and contain the node
0015 %   numbers of the triangulation of the micro sampling domain.
0016 %
0017 %   MACROVERTICES contains the node coordinates of the macro triangles.
0018 %
0019 %   CONSTRAINTS contains information about the micro boundary problem. The
0020 %   format depends on whether dirichlet or periodic boundary conditions are
0021 %   used. See make_contraints.m for details.
0022 %
0023 %   DELTA is the width or height of the square sampling domain.
0024 %
0025 %   BCTYPE is either "dirichlet" or "periodic" and determines the boundary
0026 %   condition of the micro problem.
0027 %
0028 %   A_MACRO has dimension 4 x 4. The vertices are numbered
0029 %   anti-clockwise.
0030 %
0031 %   This function should not be modified.
0032 %
0033 %
0034 %   The code is available at http://anmc.epfl.ch/ and described in
0035 %   further detail in
0036 %
0037 %   A. Abdulle and A. Nonnenmacher
0038 %   "A short and versatile finite element multiscale code for
0039 %   homogenization problems"
0040 %   Computer Methods in Applied Mechanics and Engineering,
0041 %   http://dx.doi.org/10.1016/j.cma.2009.03.019
0042 %
0043 %   Please cite this article in any publication describing research
0044 %   performed using the software.
0045 %
0046 %
0047 %   Email           : assyr.abdulle@epfl.ch and achim.nonnenmacher@epfl.ch
0048 %   Last updated    : 04/29/2009 with MATLAB 7.4
0049 %
0050 %   FE_HMM2D is Copyright (C) 2009 A. Abdulle and A. Nonnenmacher.
0051 %   The software is provided free for non-commercial use unter the terms of
0052 %   the GNU General Public License. See "copyright.m" for full details.
0053 
0054 %
0055 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0056 
0057 
0058 % generate micro-coordinates
0059 bary=sum(MacroVertices(:,:))/3;
0060 [MicroCoordinates]=micromesh_coords(NMicro, bary, delta);
0061 
0062 % initialization of variables
0063 NoOfNodes = size(MicroCoordinates,1);
0064 A = sparse(NoOfNodes,NoOfNodes);
0065 
0066 % stiffness matrix assembly
0067 % For meshes with many DOF this can be very slow due to Matlab's way of
0068 % handling sparse matrices. An optimized code can be found at the very
0069 % bottom of this file.
0070 for j=1:size(MicroElements,1)
0071     A(MicroElements(j,:), MicroElements(j,:)) = ...
0072         A(MicroElements(j,:), MicroElements(j,:)) + ...
0073         hmm_microstima_quad(MicroCoordinates(MicroElements(j,:),:),...
0074             bary, epsilon);
0075 end
0076 
0077         
0078         
0079 
0080 % solution of cell problems dilatation for the macro element
0081 D=[MacroVertices(2,:)-MacroVertices(1,:);MacroVertices(3,:)-MacroVertices(1,:)]';
0082 
0083 % get constraint matrix and corresponding rhs for the lagrange multiplier
0084 
0085 switch(lower(bctype))
0086     case{'dirichlet'}
0087         [ConstraintMat, ConstraintRhs]=...
0088             micro_constraints_dirichlet(Constraints, D, MacroVertices ,...
0089             MicroElements, MicroCoordinates, 0);
0090     case{'periodic'}
0091         [ConstraintMat, ConstraintRhs]=...
0092             micro_constraints_periodic(Constraints, D, MacroVertices ,...
0093             MicroElements, MicroCoordinates, 0);
0094 end
0095 
0096 
0097 % assemble matrix for constraint system and corresponding rhs
0098 AConstr=[A ConstraintMat';...
0099     ConstraintMat sparse(size(ConstraintMat,1),size(ConstraintMat,1))];
0100 RhsConstr=[zeros(NoOfNodes,3);ConstraintRhs];
0101 
0102 % finally solve
0103 x=AConstr\RhsConstr;
0104 alpha=x(1:size(MicroCoordinates,1),:);
0105 
0106 % area of macro-element and micro-domain
0107 K_macro=.5*det([1,1,1; MacroVertices']);
0108 K_micro=delta^2;
0109 
0110 % the entries for the macro stiffness matrix
0111 A_macro=K_macro/K_micro *alpha'*A*alpha;
0112 
0113 end
0114 
0115 
0116 
0117 
0118 %
0119 % Appendix
0120 %
0121 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0122 
0123 % % stiffness matrix assembly - fast version
0124 % indexx=zeros(4*4*size(MicroElements,1),1);
0125 % indexy=zeros(4*4*size(MicroElements,1),1);
0126 % value= zeros(4*4*size(MicroElements,1),1);
0127 %
0128 % for j=1:size(MicroElements,1)
0129 %     M=hmm_microstima_quad(MicroCoordinates(MicroElements(j,:),:),...
0130 %             bary, epsilon);
0131 %
0132 %     tmpa=repmat(MicroElements(j,:),1,4);
0133 %     tmpb=repmat(MicroElements(j,:),4,1);
0134 %     tmpM(:,:)=M;
0135 %     indexx(1+(j-1)*16:j*16)=tmpa(:);
0136 %     indexy(1+(j-1)*16:j*16)=tmpb(:);
0137 %     value(1+(j-1)*16:j*16)=tmpM(:);
0138 % end
0139 % A=sparse(indexx,indexy,value,...
0140 %             size(MicroCoordinates,1),size(MicroCoordinates,1));

Generated on Tue 21-Jul-2009 10:55:32 by m2html © 2003