• 小拓扑
    了解作者
  • matlab
    开发工具
  • 2KB
    文件大小
  • zip
    文件格式
  • 0
    收藏次数
  • 1 积分
    下载积分
  • 0
    下载次数
  • 2021-04-19 21:24
    上传日期
热传导拓扑优化得以个算例,关键在约束条件跟加热点。
topt2.zip
  • topt2.m
    5.3KB
内容介绍
% a 99 line topology optimization code by Ole Sigmund,October 199 %整块板均匀加热,四周边界给定温度T=0。 function topt2(nelx,nely,volfrac,penal,rmin) % initialize x(1:nely,1:nelx)=volfrac; % x是设计变量 loop=0; % 存放迭代次数的变量 change=1.; % 每次迭代,x的改变值,用来判断何时收敛 % start ineration while change>0.01 %当两次连续目标函数迭代的差<=0.01时,且迭代次数小于50次,迭代结束 loop=loop+1; xold=x; %把前一次的设计变量付给xold % FE analysis [U]=FE(nelx,nely,x,penal); %有限元分析,得到位移矢量U % objective function and sensitivity analysis 目标函数和灵敏度分析 [KE]=lk; % 单位刚度矩阵 c=0.; % 用来存放目标函数的变量柔度,使得刚度最大,柔度最小 for ely=1:nely for elx=1:nelx n1=(nely+1)*(elx-1)+ely; %左上角的单元节点 n2=(nely+1)*elx+ely; %右上角的单元节点 Ue=U([n1;n2;n2+1;n1+1],1); %所示单元的自由度,左上,右上,右下,左下 c=c+(0.001+0.999*penal*x(ely,elx)^penal)*Ue'*KE*Ue; % 计算目标函数的值(柔度) dc(ely,elx)=-0.999*penal*x(ely,elx)^(penal-1)*Ue'*KE*Ue; % 目标函数的灵敏度(导数) end end % filtering of sensitivities 灵敏度过滤 [dc]=check(nelx,nely,rmin,x,dc); %灵敏度过滤,为了边界光顺一点 % design update by the optimality criteria method 最佳优化方式设计更新 [x]=oc(nelx,nely,x,volfrac,dc); % print result 屏幕上显示迭代信息 change=max(max(x-xold)); % 计算x的改变量 abs绝对值,两个max指x行y列。 disp(['It.:' sprintf( '%4i',loop) ' Obj.:' sprintf(' %10.4f',c) ... %It 迭代,obj目标函数(柔度), ' Vol.:' sprintf('%6.3f',sum(sum(x))/(nelx*nely)) ... %vol材料体积比 ' ch.:' sprintf('%6.3f',change)]) % plot densities 优化结果的图形显示 colormap(gray);imagesc(-x);axis equal;axis tight; axis off;pause(1e-6); % axis坐标轴控制函数 axis off 坐标系设为不可见,(1e-6)10的-6次方 end % imagesc 函数显示灰度图像 下面的代码是具有两个输入参数的 imagesc 函数显示一副灰度图像 imagesc(1,[0,1]); colormap(gray); imagesc 函数中的第二个参数确定灰度范围。 end % FE analysis 有限元求解子程序 function [U]=FE(nelx,nely,x,penal) %自定义函数,最后返回[U] [KE]=lk; % sparse 把一个全矩阵转化为一个稀疏矩阵,只存储每一个非零元素的三个值:元素值,元素的行号和列号 % 总体刚度矩阵的稀疏矩阵 % *2是因为x,y都有一个数 K=sparse((nelx+1)*(nely+1),(nelx+1)*(nely+1)); %k表示整体刚度矩阵 F=sparse((nely+1)*(nelx+1),1); %F表示力的向量 U=sparse((nely+1)*(nelx+1),1); %U表示整体变形 for elx=1:nelx for ely=1:nely %一列列的排序,y*x n1=(nely+1)*(elx-1)+ely; %左上 n2=(nely+1)*elx +ely; %右上 % 左上,右上,右下,左下 自由度 % 一个点有两个,所以要*2。第一个从1开始,所以*2之后要-1。 edof=[n1;n2;n2+1;n1+1]; K(edof,edof)=K(edof,edof)+(0.001+0.999*x(ely,elx)^penal)*KE; % 将单元刚度矩阵组装成总的刚度矩阵 end end % define loads and supports 定义荷载和支撑(半MBB梁) F(:,1)=0.01; %施加载荷 垂直单元力 %施加位移约束 a=union([1:1:(nely+1)],[1:(nely+1):(nely+1)*nelx+1]); b=union([(nely+1)*nelx+1:1:(nely+1)*(nelx+1)],[(nely+1):(nely+1):(nely+1)*(nelx+1)]); fixeddofs =[a,b]; %固定节点 alldofs =[1:(nely+1)*(nelx+1)]; %固定所有节点 % 因无约束自由度与固定自由度的不同来找到无约束自由度 freedofs =setdiff(alldofs,fixeddofs); % solving 求解线性方程组,得到节点自由度的位移值 U(freedofs,:)=K(freedofs,freedofs)\F(freedofs,:); U(fixeddofs,:)=0; %受约束的节点自由度的位移值设为0 end % mesh-independency filter 独立网格过滤 敏度过滤技术子程序 function [dcn]=check(nelx,nely,rmin,x,dc) %dcn 清零,dcn用来保存更新的目标函数灵敏度 dcn=zeros(nely,nelx); %遍历所有单元 for i=1:nelx for j=1:nely sum=0.0; %遍历于这个单元相邻的单元 for k=max(i-floor(rmin),1):min(i+floor(rmin),nelx) %floor()函数时取整的意思 for l=max(j-floor(rmin),1):min(j+floor(rmin),nely) fac=rmin-sqrt((i-k)^2+(j-l)^2); %卷积算子(重量系数)Hf=过滤半径-控制函数(两元素中心距离) sum=sum+max(0,fac); dcn(j,i)=dcn(j,i)+max(0,fac)*x(l,k)*dc(l,k); end end dcn(j,i)=dcn(j,i)/(x(j,i)*sum); end end end % Element stiffness matrix 单元刚度矩阵 function [KE]=lk KE=[2/3 -1/6 -1/3 -1/6 -1/6 2/3 -1/6 -1/3 -1/3 -1/6 2/3 -1/6 -1/6 -1/3 -1/6 2/3]; end % optimality criteria update 优化迭代模块 function [xnew]=oc(nelx,nely,x,volfrac,dc) l1=0; l2=100000; %用于体积约束的拉格朗日乘子 move=0.2; while (l2-l1>1e-4) lmid=0.5*(l2+l1); xnew =max(0.001,max(x-move,min(1.,min(x+move,x.*sqrt(-dc./lmid))))); if sum(sum(xnew))-volfrac*nelx*nely>0; l1=lmid; else l2=lmid; end end end
评论
    相关推荐
    • HyperWorks14.0自学教材.zip
      针对Hypermesh 14.0软件的有限元前处理功能,以精选的案例为主线,介绍Hypermesh 的有限元建模和后分析软件接口的全过程,重点介绍Hypermesh 14.0各个模型功能及操作步骤,并结合案例介绍Hypermesh与主流CAE软件联合...
    • ANSYS_12[1][1].0_Workbench-分析.rar
      Ansys workbench12分析,ANSYS入门学习资料
    • 88行matlab拓扑优化代码-timoshenko-fem:提莫申科梁的有限元分析
      88行matlab拓扑优化代码有限元方法:Timoshenko梁和泊松方程的例子 该项目是编码有限元方法的一个小例子。 第一个编码示例是确定Timoshenko悬臂的模态频率。 第二种方法解决了对许多物理现象进行建模的泊松方程。 ...
    • top88r.zip
      88行热传导程序,用于散热结构的拓扑优化研究,希望帮到大家
    • toph.zip
      用于散热结构拓扑优化的研究,基于稳态热传导展开
    • HyperWorks14.0自学教材.rar
      针对Hyper mesh 14.0软件的有限元前处理功能,以精选的案例为主线,介绍Hyper mesh 的有限元建模和后分析软件接口的全过程,重点介绍Hyper mesh 14.0各个模型功能及操作步骤,并结合案例介绍Hyper mesh与主流C A E...
    • 5.rar
      针对一个两相混合物问题,定义一个用户滑动速度的UDF
    • 4.rar
      关于演示DEFINE_CHEM_STEP 的UDF的例子
    • Desktop3.rar
      用户使用fluent自定义计算气浊率的udf
    • SIM800C_MQTT.rar
      使用SIM800C模块,使用MQTT协议,连接中国移动onenet平台,能实现数据的订阅、发布、存储等