关于水平集方法的拓扑优化代码-TOPLSM.m

function TOPLSM（DomainWidth DomainHight EleNumPerRow EleNumPerCol LV LCurFEAInterval PlotInterval TotalItNum）
 %=================================================================%
% TOPLSM a 199-line Matlab program is developed and presented here for the
% mean compliance optimization of structures in 2D with the classical level set method. 
%
% Developed by: Michael Yu WANG Shikui CHEN and Qi XIA
% First Version : July 10 2004
% Second Version: September 27 2004
% Last Modification:October 31 2005 optimize the code.
%
% The code can be downloaded from the webpage:
% http://www2.acae.cuhk.edu.hk/~cmdl/download.htm
%
% Department of Automation and Computer-Aided Engineering 
% The Chinese University of Hong Kong
% Email: yuwang@acae.cuhk.edu.hk
%
%Main references:
% （1.）M.Y. Wang X. M. Wang and D. M. GuoA level set method for structural topology optimization
% Computer Methods in Applied Mechanics and Engineering 192（1-2） 227-246 January 2003
%
%（2.） M. Y. Wang and X. M. Wang PDE-driven level sets shape sensitivity and curvature flow for structural topology optimization
% CMES: Computer Modeling in Engineering & Sciences 6（4） 373-395 October 2004.
%
%（3.） G. Allaire F. Jouve A.-M. Toader Structural optimization using sensitivity analysis and a level-set method  
% J. Comp. Phys. Vol 194/1 pp.363-393 2004.  

%Parameters:
% DomainWidth : the width of the design domain;
% DomainHight : the hight of the design domain;
% EleNumPerRow : the number of finite elements in horizontal direction;
% EleNumPerCol : the number of finite elements in vertical direction;
% LV : Lagrange multiplier for volume constraint;
% LCur : Lagrange multiplier for perimeter constraint whose shape sensitivity is curvature;
% FEAInterval : parameters to specify the frequency of finite element
% analysis;
% PlotInterval : parameters to specify the frequency of plotting;
% TotalItNum : total iteration number.
%=================================================================%

% Step 1: Data initialization
EW = DomainWidth / EleNumPerRow;       %  The width of each finite element.
EH = DomainHight / EleNumPerCol;          % The hight of each finite element.
M = [ EleNumPerCol + 1  EleNumPerRow + 1 ]; % the number of nodes in each dimension
[ x  y ] = meshgrid（ EW * [ -0.5 : EleNumPerRow + 0.5 ]  EH * [ -0.5 : EleNumPerCol + 0.5 ]）; 
[ FENd.x FENd.y FirstNodePerCol ] = MakeNodes（EleNumPerRowEleNumPerColEWEH）; % get the coordinates of the finite element nodes
Ele.NodesID = MakeElements（ EleNumPerRow EleNumPerCol FirstNodePerCol ）; 
LSgrid.x = x（:）; LSgrid.y = y（:）;                           % The coordinates of each Level Set grid
 for i = 1 : length（Ele.NodesID（:1））
    Ele.LSgridID（i） = find（（LSgrid.x - FENd.x（Ele.NodesID（i1）） - EW/2）.^2 +... % get the ID of the level set grid that lies in the middle of a finite element
       （LSgrid.y - FENd.y（Ele.NodesID（i1）） - EH/2）.^2 <= 100*ep