算例(反分析).doc

同济曙光GeoFBA2D 反分析模块洞室开挖反分析问题描述：地下洞室的开挖，开挖边界与开挖影响区域取为如图所示。在开挖后测点P1（198.90,119.0）,P2(183.4,75.0),P3(214.4,103.5),开挖后通过监测发现测线1收缩0.05m，测线2收缩0.05m，测线3收缩0.08m，所用的地层参数为E200000kPa，0.4，kN/m3,c=12kPa，5kPa,k=0.46。通过量测信息，反演地层的力学参数弹性模量E。P2(183.4,103.5)（183.4，75）测线3P3(214.4,103.5)(198.90,119.0)P1（110.04，19）（301.70，181.6）测线1测线2(214.4,75)解题步骤：A. 定义材料1 从主菜单上选择：有限元分析材料；2 单击添加，将新添加的材料改名为地层；3 输入地层参数为E200000kPa，0.4，kN/m3,c=12kPa，5kPa,k=0.46。B. 绘制几何模型命令：rectangle第一点：110.04，19第二点：301.70，181.6命令：line第一点：183.4，75下一点：214.5，75下一点：214.4，103.5下一点：命令：line第一点：183.4，75下一点：183.4，103.5下一点：命令：point指定点的位置：198.90,119.0指定点的位置：183.4，103.5指定点的位置：214.4，103.5指定点的位置：从绘制菜单上选择圆弧，3点画圆弧模式。命令：arc指定圆弧的起点，或者圆心Center:单击P3点指定圆弧的第二点，或者圆心Center:单击P1点指定圆弧的终点：单击P2点。这样图形建模完毕。C. 保存文件。D. 寻找封闭面，赋土层材料单击命令 网格自动寻找封闭面，完成封闭面的寻找。选中封闭面，单击属性工具，弹出修改属性对话框。将类型改为地层。E. 施加边界条件单击工具条，施加边界条件，选中两竖直边施加水平方向上的固定约束。选择底边施加垂直方向上的固定约束。F. 定义施工步与开挖过程。单击工具条，弹出对话框定义为一个开挖步。到第1施工步，将洞室挖去。G. 定义测线定义测线前，必须先画好准备定义测线的普通直线。为此在点P1，P2，P3之间画直线。P1P3P2将施工步调到施工步1单击菜单 反分析测线。命令：mealine选择作为测线的直线或圆弧线段，及作为测点的点：单击P2,P1,Line1,删除原线段？Yes/NoY在选择时也可在工程窗口中选中对应的点与直线通过同样的步骤可以定义其他的两条测线。完成测线定义后的工程窗口测线2测线1测线3修改测线1的属性，选中测线1，单击工具条，弹出属性对话框，单击测量组数，详细信息，弹出如下对话框。选中施工步1添加到当前测量组。单击测量数据，详细信息，弹出如下对话框。在测点1处输入量测信息为收缩0.05m。选中测线2与测线3重复上述步骤，修改测线信息。测线2收缩0.05m，测线3收缩0.08m。保存文件。H. 定义反分析起始步与终止步单击 反分析反分析起始步与终止步。I. 反分析控制定义反分析目标为地层的弹性模量E。J. 反分析求解单击网格全自动生成网格。单击反分析反分析求解。出现反分析求解过程跟踪图， 反分析跟踪图与所选的优化分析方法密切相关。跟踪图会因所选的分析方法而不同。K. 查看计算结果计算完毕后在命令窗口，可以看到反分析的结果。 原始值 反演值由反演值计算的预测值实测值8 反分析计算数据文件格式说明反分析计算过程中，除生成.DAT文件外（文件格式见同济曙光有限元正分析手册），还将生成.MEA和.BAA文件，如果反分析优化方法不是单纯形法，还有.INP文件。以下对这些文件格式进行说明。声明：计算文件格式仅供用户参考，由于计算文件是由同济曙光软件自动生成，因此实际使用的计算文件格式可能与以下说明的格式有所不同！（1）.BAA文件数据输入格式示例BACK(analysis)*READ(11,(A1,A)LL,TXT* * READ(TXT,*) K,(IDR(K,I),I=1,4),VAR0(K,1),VAR0(K,2)* * IF (IBFLAG.GT.0) THEN ！INPUT BOUNDARY LOADS PARAMETERS * READ(11,(A1,A)LL,TXT* READ(TXT,*) K,(NBCXY(K,I),I=1,5),BCXY(K,J),J=1,6)* ENDIF*-+-+-+-+-+-+-+-+-+ 1 1 2 1 0 60000.00 0.0 2 1 3 1 1 5500.00 1500.0E 3 1 5 1 1 55900.00 30000.0相应变量名称：K,(IDR(K,I),I=1,4),VAR0(K,1),VAR0(K,2)K,(NBCXY(K,I),I=1,5),BCXY(K,J),J=1,6)变量含义：LL输入结束标志（=E或e 表示输入结束，否则继续输入）K输入行编号构造地应力反演参数初始值，与CXY(I)中各值一一对应。IDR(K, 1) 反演种类项标志 （=0 反演地应力CXY(J)参数 =1 反演材料AMT(I,J)内参数； =2 反演BCXY(I,J)内边界荷载分布参数）IDR(K,2)为各参数项反演标志 （=1 反演； =0 已知）IDR(K, 3)、IDR(K, 4)需反演参数在CXY(J)、AMT(I,J)、BCXY(I,J)数组内的行I和列J，地应力反演中，IDR(*, 3)均填1（行），IDR(*, 4)值（列）对应cxy(*)中的元素VAR0(K,1)为反演参数初始值（IMTHD=1，2），或反演参数平均值（IMTHD=3，4，5）VAR0(K,2)为反演参数初始值范围（IMTHD=1，2任意值，IMTHD=3，4，5需要）以下反演边界荷载参数（存在IDR(K, 1)=1且IDR(K, 4)=1）时输入：NBCXY(K,1)、NBCXY(K,2) 边界荷载存在之开挖步、增量步NBCXY(K,3)、NBCXY(K,4) 边界荷载两端（分布荷载，中间可有结点）结点号NBCXY(K,5) 边界荷载方向（=1 总体坐标，=2 下图所示局部坐标）BCXY(K,6) 分别为当前步首结点沿荷载线的法向、切向、转动向和末结点沿荷载线的法向、切向和转动向荷载（kN）结点2结点1注：第一列以 * 开始的均为说明（2）.MEA文件数据输入格式示例：第一组数据MEASURE DATA1 2 0.0 0.0 0 1 1 1 相应变量名称 MEASTEP, IFLAG, X0, Y0, ISTAG1,ISTEP1,ISTAG2,ISTEP2变量含义： MEASTEP量测次数 IFLAG量测点后处理显示标志，一般赋2 X0,Y0量测点相对坐标系原点的坐标 ISTAG1,ISTEP1有量测值的起始开挖步与相应增量步，一般ISTEP1=本次开挖最大增量步（注意：ISTAG1,ISTEP1一般填（0，1），ISTAG1=0时，以下该步量测信息可不必填，因该初始步位移、内力、应变计算值为0，但初始应力量测信息计算与实测均非0，应填） ISTAG2,ISTEP2有量测值的终止开挖步与相应增量步，一般ISTEP2=本次开挖最大增量步注：采用量测值增量进行反演，即ISTAG2,ISTEP2对应量测值减去ISTAG1,ISTEP1对应量侧值。第二组数据 共MEASTEP组数据输入格式示例： 10 1 1 0 5 2 2 2 2 2 0 0 0 0相应变量名称NO,ISTAGE,ISTEP,NT1,(MPOI(I,1),I=1,NT1),NT2,(MPOI(I,2),I=1,NT2),NT3,NT4,NT5, NT6变量含义：NO量测点数ISTAGE,ISTEP当前量测数据对应的开挖步、增量步NT1,(MPOI(I,1),I=1,NT1)绝对位移测线数,各测线对应的测点数NT2,(MPOI(I,2),I=1,NT2)相对位移测线数（相对点为各自测线中的最大号测点），各测线对应的测点数NT3梁单元弯矩测点数NT4结构轴力测点数NT5应变测点数（实体元和接触元）NT6应力测点数（实体元和接触元）注意：弯矩、轴力、应变和应力按测点计，无测线（或按测线理解可认为每个测线仅一个测点）第三组数据 共 NO 组 数据输入格式示例：1 19.821477 96.280000 1 0.006210 0.000000 10 12 24.063189 96.280000 1 0.000000 0.000000 10 13 19.981500 94.780000 1 0.006453 0.000000 10 14 24.092717 94.780000 1 0.000000 0.000000 10 15 21.990875 98.518870 1 0.000000 45.903000 10 16 19.821477 96.280000 1 0.004887 45.903000 10 17 21.990875 98.518870 1 0.000000 61.745000 10 18 19.981500 94.780000 1 0.005265 61.745000 10 19 21.990875 98.518870 1 0.000000 119.3422500 10 110 24.092717 94.780000 1 0.003105 119.342250 10 1相应变量名称 K1, XXI, YYI, K2, SDIS, ALFA1, K3, K4变量含义：K1测点编号XXI,YYI测点相对坐标K2坐标标志， =1 X,Y 坐标 =2 r,q 极坐标SDIS绝对位移、相对位移、弯矩、轴力、应变或应力ALFA1位移、应变或应力方向。位移：360转角位移；对弯矩、轴力可取任意值；对应变：=1，2，3 应变在x, y (接触面切向和法向)和xy向；对应力：=1，2，3，4 应力在x, y (接触面切向和法向)、xy向和z向。K3量测点所在单元材料性质类型（与IAMT(I)相对应）K4量测值性质标志， =1 表示量测位移为累计量 =2 表示量测位移为增量（对开挖步增量步而言）注：关于相对位移SDIS和测线方向角ALFA1填写的说明：相对位移测线中，每根测线最下面测点为基点，其他测点测值投影到测线上再减去基点在测线上的投影值（同一测线角度相同），当然还得考虑施工步和增量步的差值。即该投影差值在前后两不同施工步和增量步有两个值，这两个值的差值为反演时采用量测值，其值与计算值比较，求目标函数。单位与网格单位一致。相对位移SDIS（设i点记为di）一般填收敛值，设1根相对位移测线上共有1.in个测点，i为任一测点，n点为基点，位移收敛值D(di-dn)填正值，测线方向角ALFA1为测线矢量in与x轴夹角；位移收敛值D(di-dn)表示拉伸时填负值，测线方向角ALFA1同上。（3）.INP文件注1：当优化方法为单纯形法时，无此参数文件。注2：以*号开头的行为注释行。遗传算法：* heredity algorithm option40 12 0.02 120 0.5 1 1 16 1540 为种群规模；12 需要反演的参数个数；0.02变异概率；120迭代的最大代数；0.5杂交概率；1为均匀杂交（0为单点杂交）；1为双亲一子（0为双亲两子）；16为迭代控制精度（216）；15为迭代终止条件（连续15迭代最优适应值没有进化终止）。模拟退火算法：* simulate annealing algorithm option10000 0.9810000为初始温度；0.98降温系数。阻尼最小二乘法：* damping least square algorithm option100 5000 0.0001100为最大迭代代数；5000为控制条件（目标函数小于5000，则终止）；0.0001如所有参数的增值绝对值的和小于所有参数绝对值的和