【COMSOL】全流程实践与解析 Demo6-THMD热-流-固-损伤耦合

发布时间:2026/9/1 6:26:23
【COMSOL】全流程实践与解析 Demo6-THMD热-流-固-损伤耦合 目录一、模型向导二、建模1. 全局参数2. 定义1插值2变量3状态变量4分段函数3. 几何4. 材料5. 物理场1固体力学2达西定律3多孔介质传热6. 网格三、计算1. 求解器设置2. 计算损伤及属性效果一、模型向导空间【二维】物理场【多孔弹性 实体】【多孔介质传热】研究【稳态】向导完成后再添加一个【瞬态】二、建模1. 全局参数名称表达式值描述fc160[MPa]1.6E8 Pa单轴抗压强度phi40[deg]0.69813 rad内摩擦角poro_00.010.01初始孔隙率poro_r1.00E-030.001残余孔隙率mu0.20.2泊松比k01e-19[m^2]1E-19 m²初始渗透率pinj30[MPa]3E7 Pa注入压力定值landa04[W/m/K]4 W/(m·K)初始固体热传导系数Cs0950[J/kg/K]950 J/(kg·K)固体比热容h_tc01000[W/(m^2*K)]1000 W/(m²·K)初始传热系数Tinit373.15[K]373.15 K初始温度Tinj293.15[K]293.15 K注入温度h_expan_s6e-6[1/K]6E-6 1/K固体热膨胀系数h_expan_l1e-3[1/K]0.001 1/K用于算c2h_expan_m6.2e-6[1/K]6.2E-6 1/K用于算c22. 定义1插值【定义】→【函数】→【插值】创建两个插值函数一个是初始杨氏模量 E另一个是初始抗拉强度 ft【数据源】文件选择 MATLAB 随机生成的 E 和 ft 数据表后导入MATLAB 参考代码适用于2019以下版本value 16; %抗拉强度 ft 参考数值 %value 35; %杨氏模量 E 参考数值 hang 0.3; %m %宽度 lie 0.3; %m %高度 n 0.002; %取样间隔 pro wblrnd(value,20,hang/n,lie/n); mat []; A linspace(0,hang,hang/n); %矩阵转换 mat ones(hang/n,1)*A; B mat(:); C mat; C C(:); D pro; D D(:); T [B C D]; min(D) max(D) mean(D) xlswrite(weibull-ft.xlsx,T); %xlswrite(weibull-E.xlsx,T); [x1,y1]meshgrid(A,A); surf(x1,y1,pro) shading interp杨氏模量 E 的单位为 GPa抗拉强度 ft 单位为 MPa2变量名称表达式F1solid.sp1-ft0F2-solid.sp3(1sin(phi))/(1-sin(phi))*solid.sp1-fcE0int1(x,y)et0ft0/E0*0.5ec0fc/E0poroporo_r(poro_0-poro_r)*exp(alpha*p)alpha3*(1-2*mu)/E0/(poro_0-poro_r)kmk0*(poro/poro_0)^3*exp(5*D)Dtif((F10)(d(F1,TIME)0),min(0.99,max(0,1-(et0/solid.ep1)^2)),0)Dc-if((F20)(d(F2,TIME)0),min(0.99,max(0,1-(ec0/solid.ep3)^2)),0)DDmin(0.99,Dtabs(Dc))ft0int2(x,y)landalanda0*exp(D*5)CsCs0h_tcif(abs(D)0.95,5*h_tc0,h_tc0)biotif(D0.9,0.5,1)c1biotc2poro*h_expan_l(1-poro)*h_expan_s-h_expan_m*(1-biot)Qssolid.K*h_expan_s*Td*d(solid.eth11solid.eth22,TIME)Tdsolid.Tdiff3状态变量先在【显示更多选项】中勾选上【方程贡献】comsol 6.4版本是勾选这个其他版本可能不是在解释里找【定义】→【方程贡献】→【状态变量】表达式if(DDD,DD,D)4分段函数【定义】→【函数】→【分段】3. 几何【正方形】边长 0.3【圆】半径 0.015xy 0.15【布尔】→【差集】正方形减去圆形【形成联合体】4. 材料1空材料缺失属性会在设置完物理场的边界条件后自动显示可以先创建好边界条件后再来设置具体的属性值2water添加【比热率】属性这里设置为 1.00665. 物理场1固体力学【辊支承】左 下边界约束法向位移【边界荷载】上 右 内部圆3个分开竖向荷载小右边界内部边界【线弹性材料】→【热膨胀】体积参考温度 Tref 设置为 Tinit温度用 ht2达西定律【流体】mat2 指的是材料中的水【初始值】【压力】内部边界设置压力设置为 0.1 MPa【入口】内部边界设置压力【质量源】耦合项3多孔介质传热初始温度用 Tinit【流体】【基体】【温度】四周边界设置【初始值】和【温度】中的 T 设置为 Tinit【热通量】内部边界设置【热源】耦合项6. 网格【自由三角形】→【大小 1】整个域单元最大 0.003最小 0.0004最大增长率 1.05【大小 2】内部边界单元最大、最小0.001【细分方法】Delaunay会根据几何形状和网格大小要求自动把计算区域划分成质量较好的三角形网格尽量避免生成又细又扁的三角形三、计算1. 求解器设置【稳态】求解器禁用【固体力学】的【压力】、【达西定律】的【入口】、【多孔介质传热】的【热通量】【瞬态】求解器设置【输出时步】、【相对容差】禁用【达西定律】中的【压力】配置【因变量值】打开【显示默认求解器】【时间步进】设置【求解时显示】设置【全耦合】设置2. 计算计算【稳态】【应力】改为【损伤】修改表达式禁用【变形】删了也行【温度】的单位改成摄氏度 degC【损伤】、【压力】、【温度】的数据集改为研究2计算【瞬态】计算过程中会展示损伤变化【损伤】计算结果【压力】【温度】【孔隙率】