matlab建模与仿真实例_MATLAB辅助ANSYS APDL进行负泊松比铰链结构建模仿真

       ANSYS有着强大的参数化建模能力，如果对命令流不太熟悉，或者需要建模复杂，用命令流写代码计算太麻烦，有时候需要借助matlab强大的矩阵运算能力。本例采用matlab辅助ansys进行负泊松比结构的参数化建模，如图所示，选择中间square型结构，黑色杆和灰色杆分别是模量极大的直角杆。在B点两直角杆的角点铰接， 在A点和D点端点铰接。

    图片引自：https://doi.org/10.1098/rspa.2014.0538

    为了保证数值稳定性，铰接小模量的弹性杆(蓝色)。画出示意图分别对关键点编号(黑色编号)角点铰接处有两个关键点分别属于两个直角杆、线编号(品红色)。

    很容易找到规律。并以左下角单元为基本单元，进行排列展开。在matlab中求解关键点的位置坐标和连线关系矩阵。进行验证后，便完成了任意行任意列的阵列排列命令流建模。

       在anays中进行梁单元的铰接模拟需要用到以下命令流。有两种，分别是自由度耦合命令cp与自由度释放命令endrelease，关于这两个命令的简单案例，可以查阅王新敏的《ANSYS工程结构数值分析》p350~p353。根据该问题的特点，选择在直角点铰接处使用endrelease命令，在三杆连接处使用cp命令。

    在matlab中运行结尾的程序，便可生成ansys可以直接运行的命令流文件array.txt，命令流文件中提取了运算结果的位移数据等，以便在matlab中计算其等效泊松比。具体公式为，见论文：

https://doi.org/10.1098/rspa.2014.0538

    生成命令流文件array.txt后，复制到ansys的当前工作目录下(程序最后注释的部分可以将命令流文件自动发送到指定的工作目录下)，在命令流输入窗口直接输入array.txt回车即可，完成该负泊松比结构的仿真计算。

     软件安装目录和软件安装管家最近都无法安装软件，如果有安装matlab的需求，可以选择这个网站，里面有安装包和详细的安装教程：

    http://www.softgj.com/win/matlab

    有的版本安装包装上后，会出现一些问题。之间试过2018的一个安装包，打开后初始化特别特别慢。最后换了2019a就好了，所以建议安装合适的版本。

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%求解任意m行任意n列的square布周期性结构排布周期性结构的ansys建模相关设置的命令流%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%clear allclcclose all%% square % m行 n列m = 2;n = 3;%% gamma ga = 70; %角度制 gamma = ga*pi/180;%弧度制theta = 90*pi/180;alpha = 90*pi/180;p = 1;%mmOB = p*sin(gamma)/sin(alpha/2);%% 坐标coor_release(1,:) = [OB*cos(alpha /2)-p*cos(gamma),0];coor_release(2,:) = [OB*cos(alpha /2) + p*cos(gamma),0];coor_release(3,:) = [OB*cos(alpha /2) , OB*cos(alpha /2)];coor_release(4,:) = [0 , OB*cos(alpha /2) - p*cos(gamma)];coor_release(5,:) = [0 , OB*cos(alpha /2) + p*cos(gamma)];coor_release(6,:) = [-coor_release(3,1)  coor_release(3,2)];coor_release(7,:) = [-coor_release(1,1)  coor_release(1,2)];coor_release(8,:) = [-coor_release(2,1)  coor_release(2,2)];coor_release(9,:) = [-coor_release(3,1)  -coor_release(3,2)];coor_release(10,:) = [coor_release(4,1)  -coor_release(4,2)];coor_release(11,:) = [coor_release(5,1)  -coor_release(5,2)];coor_release(12,:) = [coor_release(3,1)  -coor_release(3,2)];dis_cell = 2*OB*cos(alpha /2);%% set 1  num 1:2+4ncoor_set1(1,:) = coor_release(9,:);coor_set1(2,:) = coor_release(9,:);coor_set1(3,:) = coor_release(10,:);coor_set1(4,:) = coor_release(11,:);coor_set1(5,:) = coor_release(12,:);coor_set1(6,:) = coor_release(12,:);num_set1 = 7;for i = 1:n-1    coor_set1(num_set1,:) = coor_set1(3,:)+i*[dis_cell 0];        num_set1 = num_set1 + 1;    coor_set1(num_set1,:) = coor_set1(4,:)+i*[dis_cell 0];        num_set1 = num_set1 + 1;    coor_set1(num_set1,:) = coor_set1(5,:)+i*[dis_cell 0];        num_set1 = num_set1 + 1;    coor_set1(num_set1,:) = coor_set1(6,:)+i*[dis_cell 0];        num_set1 = num_set1 + 1;end%% set 2coor_set2(1,:) = coor_release(8,:);coor_set2(2,:) = coor_release(7,:);coor_set2(3,:) = coor_release(1,:);coor_set2(4,:) = coor_release(2,:);num_set2 = 5;for i = 1:n-1    coor_set2(num_set2,:) = coor_set2(3,:)+i*[dis_cell 0];       num_set2 = num_set2 + 1;    coor_set2(num_set2,:) = coor_set2(4,:)+i*[dis_cell 0];       num_set2 = num_set2 + 1;end%% 阵列坐标num_array = 1;for i = 1:m %% set1 复制    for k = 1:(2+4*n)        coor(num_array,:) = coor_set1(k,:)+(i-1)*[0 dis_cell];        num_array = num_array + 1;    end  %% set2 复制      for k = 1:(2+2*n)        coor(num_array,:) = coor_set2(k,:)+(i-1)*[0 dis_cell];        num_array = num_array + 1;    endend %% 封顶    for k = 1:(2+4*n)        coor(num_array,:) = coor_set1(k,:)+m*[0 dis_cell];        num_array = num_array + 1;   enddisp('keypoint求解完成')%% line % m行 n列max_linehao = n*5*(m+1) + 5*m*(n+1);%% line haonum_point_begin_set1 = 1;num_point_pre = (2+4*n) + (n+1)*2;%set1 和set2num_point_pre_set1 = (2+4*n);%%set1行num_line = 1;for i = 1:m      num_point_begin = 1 + (i-1)*num_point_pre;      k = num_point_begin;%      %% line_hao_set1       line_hao_set1(1,:) = [k k+3];       line_hao_set1(2,:) = [k+3 k+5];       line_hao_set1(3,:) = [k+2 k+3];       line_hao_set1(4,:) = [k+1 k+2];       line_hao_set1(5,:) = [k+2 k+4];     %% set 1       for j = 1:n            for index = 1:5             line_hao(num_line,:) = line_hao_set1(index,:)+(j-1)*[4 4];%杆1             num_line = num_line + 1;           end       end     %% set 2       for j = 1:n+1         point_1 = k+num_point_pre + 4*(j-1);        point_2 = k+num_point_pre+1 + 4*(j-1);        point_3 = k+num_point_pre_set1+1 + 2*(j-1);        point_4 = k+num_point_pre_set1 + 2*(j-1);        point_5 = k + 4*(j-1);        point_6 = k+1 + 4*(j-1);           line_hao(num_line,:) = [point_1 point_4];%杆1             num_line = num_line + 1;           line_hao(num_line,:) = [point_4 point_6];%杆2             num_line = num_line + 1;          line_hao(num_line,:) = [point_3 point_4];%杆3             num_line = num_line + 1;           line_hao(num_line,:) = [point_2 point_3];%杆4             num_line = num_line + 1;           line_hao(num_line,:) = [point_3 point_5];%杆5             num_line = num_line + 1;       endend%% 封顶  %% line_hao_set1      num_point_begin = 1 + m*num_point_pre;      k = num_point_begin;%      %% line_hao_set1       line_hao_set1(1,:) = [k k+3];       line_hao_set1(2,:) = [k+3 k+5];       line_hao_set1(3,:) = [k+2 k+3];       line_hao_set1(4,:) = [k+1 k+2];       line_hao_set1(5,:) = [k+2 k+4];     %% set 1       for j = 1:n            for index = 1:5             line_hao(num_line,:) = line_hao_set1(index,:)+(j-1)*[4 4];%杆1             num_line = num_line + 1;           end       endhao = 1:size(line_hao,1);num_line_hao = [hao' line_hao];disp('line求解完成')%% 画图验证    disp('画图');    plot(coor(:,1),coor(:,2),'b*');    axis equal     hold on     %% point number     for i = 1:size(coor,1)         hao(i,1) = i;         text(coor(i,1),coor(i,2),[num2str(hao(i,1))],'VerticalAlignment','bottom','FontSize',6);        %   pause(0.7)    end%        %% line number          tanhuang_hao = 3:5:(max_linehao-2);         for i = 1:size(line_hao,1)        line(coor(line_hao(i,:)',1),coor(line_hao(i,:)',2))      end         for i = 1:size(tanhuang_hao,2)        line(coor(line_hao( tanhuang_hao(i),:)',1),coor(line_hao( tanhuang_hao(i),:)',2),'Color','cyan')       end    hold off           for i = 1:size(line_hao,1)         mid_coor = [(coor(line_hao(i,1),1)+coor(line_hao(i,2),1))/2  (coor(line_hao(i,1),2)+coor(line_hao(i,2),2))/2];         hao(i,1) = i;           text(mid_coor(1,1),mid_coor(1,2),[num2str(hao(i,1))],'Color','magenta','FontSize',8);     end  hold off fid = fopen('array.txt','w');  %% 材料截面信息 s_pre = {'finish'           '  /clear'           ' /prep7'           ' !beam单元'          ' ET,1,BEAM188'          '!材料1'           'mp,ex,1,72e6'          ' mp,prxy,1,0.3'          '!材料2'          ' mp,ex,2,1'          ' mp,prxy,2,0.49'     '     !截面信息 矩形截面 B=H'         ' SECTYPE,   1, BEAM, RECT, , 0   '       '   SECOFFSET, CENT '         ' SECDATA,0.1,0.1,0,0,0,0,0,0,0,0,0,0'}; for i = 1:size( s_pre,1)       fprintf(fid, '%s\n',s_pre{i,1}); end %% 获得keypoint命令流for i = 1:size(coor,1)  s_kpoint{i,1} = strcat('k',',',num2str(i),',',num2str(coor(i,1)),',',num2str(coor(i,2)));end%% 打印keypoint命令流   fprintf(fid, '%s\n',  '!keypoint命令流');for i = 1:size(coor, 1)    fprintf(fid, '%s\n', s_kpoint{i,1});end  %% 获得Line命令流    fprintf(fid, '%s\n',  '!line命令流');for i = 1:size(line_hao,1)  s_Linehao{i,1} = strcat('l',',',num2str(line_hao(i,1)),',',num2str(line_hao(i,2)));end%% 打印Line命令流for i = 1:size(line_hao, 1)    fprintf(fid, '%s\n',  s_Linehao{i,1});end%% 打印赋材料命令流 fprintf(fid, '%s\n',  '!赋材料');fprintf(fid, '%s\n',  'allsel'); s_material_1 = strcat('lsel',',','s',',,,',num2str(1),',',num2str(max_linehao));fprintf(fid, '%s\n',  s_material_1);fprintf(fid, '%s\n',  'latt,1,,1,,,,1'); s_material_2 = strcat('lsel',',','s',',,,',num2str(3),',',num2str(max_linehao-2),',',num2str(5));fprintf(fid, '%s\n',  s_material_2);fprintf(fid, '%s\n',  'latt,2,,1,,,,1');%% 打印mesh 命令 s_mesh = {'allsel'             'lesize,all,,,3'            'lmesh,all'}; for i = 1:size( s_mesh,1)       fprintf(fid, '%s\n',   s_mesh{i,1}); end disp('mesh完成') %% 打印release 命令 fprintf(fid, '%s\n',  '!定义release自由度约束');fprintf(fid, '%s\n',  'allsel');for k_line = 1:5:max_linehao-4   s_release = strcat('lsel',',','s',',,,',num2str(k_line),',',num2str(k_line+2));     fprintf(fid, '%s\n', s_release);      fprintf(fid, '%s\n',  'esll');       fprintf(fid, '%s\n',  'endrelease,,,rotz');       fprintf(fid, '%s\n',  'allsel');   s_release = strcat('lsel',',','s',',,,',num2str(k_line+2),',',num2str(k_line+4));     fprintf(fid, '%s\n', s_release);      fprintf(fid, '%s\n',  'esll');       fprintf(fid, '%s\n',  'endrelease,,,rotz');       fprintf(fid, '%s\n',  'allsel');end disp('release完成') %% 打印cp 命令 num_point_pre = (2+4*n) + (n+1)*2;%set1 和set2num_point_pre_set1 = (2+4*n);%%set1行 fprintf(fid, '%s\n',  '!定义cp自由度约束');fprintf(fid, '%s\n',  'allsel');for i = 1:m+1    k_line_start = 1+(i-1)*num_point_pre;    k_line_end = k_line_start + num_point_pre_set1 - 2;    for k = k_line_start:4:k_line_end         s_cp = strcat('ksel',',','s',',,,',num2str(k),',',num2str(k+1));         fprintf(fid, '%s\n', s_cp);         fprintf(fid, '%s\n',  'nslk,s');         fprintf(fid, '%s\n',  'cp,next,ux,all');         fprintf(fid, '%s\n',  'cp,next,uy,all');    endend disp('cp完成')%% 打印边界约束命令fprintf(fid, '%s\n',  'allsel');for i = 1:m    k_line_start = 1+(i-1)*num_point_pre;    k_bc_start = k_line_start + num_point_pre_set1;    k_bc_end = i*num_point_pre;    %% 左边     s_bc = strcat('dk',',',num2str(k_bc_start),',','UX');         fprintf(fid, '%s\n', s_bc);     s_bc = strcat('dk',',',num2str(k_bc_start),',','UZ');         fprintf(fid, '%s\n', s_bc);         if i==1            s_bc = strcat('dk',',',num2str(k_bc_start),',','UY');            fprintf(fid, '%s\n', s_bc);         end     %% 右边         if i==1          s_bc = strcat('dk',',',num2str(k_bc_end),',','UY');         fprintf(fid, '%s\n', s_bc);     end      s_bc = strcat('dk',',',num2str(k_bc_end),',','UX',',','0.1');         fprintf(fid, '%s\n', s_bc);        s_bc = strcat('dk',',',num2str(k_bc_end),',','UZ');         fprintf(fid, '%s\n', s_bc);         end disp('BC 完成')%% 求解后处理   s_solu_post = {'allsel'            'finish'            '/solu'            'solve'           ' /POST1 '             '/EFACET,1'             ' PLNSOL, U,SUM, 1,1.0'};               for i = 1:size(s_solu_post,1)       fprintf(fid, '%s\n',s_solu_post{i,1}); end%% 打印输出位移云图的命令流load xianshi s_output_result = xianshi;       for i = 1:size( s_output_result,1)       fprintf(fid, '%s\n', s_output_result{i,1}); end%% 打印提取节点坐标和节点位移 命令流load readresult  s_read_result =  readresult ;       for i = 1:size( s_read_result,1)       fprintf(fid, '%s\n',s_read_result{i,1}); end disp('命令流打印完成')fclose(fid);%% 复制命令流文件到ansys workspace% srcpath=[]; % 源文件全路径名% dstpath='D:\ansys\workspace\ouhe\';% 目标文件全路径名% % filename='array_cp_release.txt';% filename='array.txt';% fdt = dir([srcpath filename]);% for k = 1:length(fdt)%     fn = fdt(k).name;%     if ~isdir(fn)       % 检验是否目录，是目录跳过%         ofn = [srcpath fn];             % 源文件全路径名%         nfn = [dstpath fn];      % 目标文件全路径名%         copyfile(ofn,nfn);          % 拷贝文件%     end% end%  disp('命令流发送完成')