【Dino结构笔记】结构工程的造山_不规则壳元的建模

以下分享工程师遇到问题解决问题的过程：

现在越来越多的建筑师与结构工程师可以通过RHINO+GH对复杂曲面进行建模，通过参数化建模最后把模型导入SAP2000等软件进行结构分析。去年有幸参与一个主题乐园项目的结构分析，从中遇到了一种特别的情况，需要把不规则壳元进行有限元分析。事情是这样子的：

目前的建筑设计已经大量应用了3D打印及3D扫描技术，建筑师或设计师，只需要建立缩尺模型，采用高 精度的3D扫描技术，把模型数码化，在电脑的数码格式，可以把缩尺模型放大至实际建筑物或景观的大小，这也就是现在的人造景观的技术。目前已经可以做到“愚公造山”，结构工程师只需要设计一个钢结构  框架，如图所示，外表面可以采用GRC混凝土（带一定受拉承载力的混凝土），这种混凝土往后发展可以成为3D打印结构的主要材料。

钢结构与外蒙皮结构的关系


3D 扫描技术在建模中的应用

岩石山体扫描后数码化的效果

落在手上的问题是，需要对3D扫描的GRC混凝土壳进行有限元分析，与传统的参数化生成的壳元不同，3D扫描得到的壳元的点分布是杂乱无章的（3D扫描是得到点云，Point Cloud，这与激光的频率，点与点的距离有关，没有什么规律），需要把壳元规划划分后导入SAP2000进行计算。因为这个项目少说也有几百块壳体，批量化的划分单元是非常重要的。结构工程师遇到复杂模型的建模问题，一般通过编程解决，与传统的网格划分技术不同，这次是把网格规则化。在这个过程中，享受采用小学二年级学习的空间几何知识解决问题的快感。


以上外壳为例，在3DSMAX打开壳元的实际网格如图所示。

自编小程序用于批量划分网格 3D surface Mesh. 

1）编制程序读取3D扫描得到的不规则曲面
2）生成正交规范的点坐标
3）求解点坐标到三角面内的映谢，得到正交点的Z轴坐标。
4）最后通过正交点的三维坐标形成新的规则曲面
5）另存为DXF，导入SAP2000进行静力分析，如图所示。
6）由于有很多块壳元，写个LISP就可以自动操作起来了
7）特别三轴复杂的曲面，需要人为处理一下

规划划分后的SAP2000 模型

这里用到的小学二年级的空间知识居然是已知三点坐标求得平面方程，如下图所示。有时候简单的数学拯救了工程师。这个是对去年的有限元网格划分在工程中应用的一次小总结。现在有了Rhino与GH,这种问题又有了新的解决方法，但是对于大量单元的网格划分，如有代码进行后台计算的效率，而且还可以打包成容易入手的小工具。

 
最后发一张照片留念一下。


以上就是这篇笔记的全部内容。
以上操作用到的AutoCAD的二次开发，SAP2000的二次开发，三维网格划分的知识，空间几何的知识等。

==============================

这个是我的微信公众号[Dino结构笔记]：DINOSTRU

也可以通过扫二维码进行关注一下。

广东省地区与香港地区建筑结构风荷载规范的比较研究

小编：辛展文，指导老师：陈学伟、杨易

【1】研究背景与研究思路

风荷载设计是建筑结构设计中的一项重要内容，对沿海超高层建筑尤其重要。高层建筑抗风设计方法：规范设计、风洞试验、数值风洞模拟仿真，其中规范是基准。


2019年，香港屋宇署新修订了风荷载规范。广东省风荷载设计规范按照《建筑结构荷载规范》（2012版）实行。掌握这两地的风荷载规范，对粤港两地的工程项目非常重要。在这样的背景下，确定了本研究的选题。

研究的思路如下：

【2】两地规范对比

基本计算公式对比如下：

可以看出两地顺风向风荷载的计算形式相似，但在横风方向上，广东规范通过计算各楼层风荷载来衡量风作用，而香港规范则通过基底弯矩来衡量。计算得到横风基底弯矩后，香港规范需要根据横风基底弯矩对顺风向风荷载进行修正，广东规范则没有这部分内容。对于扭转风荷载的计算，香港规范采用一个更简单的计算模型。它计算风振加速度的公式也较为简单。

值得注意的是，广东规范只需要计算两个方向（X、Y）的风作用，而香港规范需要验算四个方向（±X1、±X2），而且香港规范可以计算任意风向

部分计算参数对比如下：

可以看出两套规范下的计算参数在取值规定与测量方法上也各有差异。计算公式和计算参数的差异性最终导致两地风荷载规范计算结果有较大的差别。

【3】风荷载计算软件编写

在掌握了规范条文后，使用面向对象编程语言DELPHI针对两地规范分别编写风载计算软件。该计算软件可以在不建模的情况下，通过输入建筑尺寸参数、结构动力特性参数、风气候特征等参数，对建筑结构进行风荷载分析。

【4】计算软件准确性检验

利用编写得到的两个计算软件，对一67层的超高层建筑进行风荷载计算。该建筑总高度为295.1m，建筑宽度为46.9m，进深为46.9m，建筑平面体型为矩形。

广东风荷载计算软件计算得到的电算结果通过PKPM进行检验，对比结果显示PKPM与自编软件计算得到的顺风向风荷载相差1.86%，横风向风荷载只相差0.34%。

香港风荷载计算软件通过手算屋面层风荷载作用的方式进行校核，最后得到的手算结果也和电算结果相符合。

【5】实例计算结果分析

对两个软件计算得到的数据进行对比，可以看出香港规范的计算结果比广东规范的要大很多。它们之间的偏差会随着层高的增加而增大，在顶层的偏差达到最大值，并且可以看出，Y方向上的偏差比X方向上的小。

为了更深入地理解这些计算结果的差异性，该研究对计算过程中涉及到的一个重要参数，不同高度处的风压值，进行计算分析，取最接近香港的城市：深圳，所采用的基本风压0.90（对应回归周期100年），为分析对象，计算得到的风压值也比香港规范的要小得多，其中最接近香港规范的是A类地貌。

【6】总结

风荷载规范是建筑结构抗风设计的重要基准，涉及建筑结构的安全和经济。本文在深入研究和详细对比粤港两地规范基础上，依据规范中的风荷载计算模型，通过面向对象语言DELPHI编程，成功开发了两款风荷载计算软件。

通过对一栋超高层建筑工程实例的风荷载分析，验证了风荷载计算软件的准确性，计算精度达到商用软件PKPM标准。

通过规范对比和工程算例发现，两地风荷载规范的计算模型和基本参数存在差异；由此导致粤港临近地区，即使建筑体型和结构动力特性相似的建筑，其风荷载规范取值存在显著差异。如本例中，顺风向风荷载基底总剪力按照广东省规范计算比香港地区的小64.98%。

香港风荷载程序的介绍与使用：

bilibili网站：https://www.bilibili.com/video/bv1Az4y1Q7uy

[Dino笔记]风荷载规范中的SFL图表问题 


副标题是：横风向广义力功率谱SFL参数的计算分析

在大部分的工程师看来风荷载的计算是一件容易的事情，就是PKPM或YJK点一个按钮的功夫，但是风荷载背后
的计算公式是相当复杂的，2012年的风荷载规范又把这个难度提高了，里面都是一些复杂的公式，如果说结构工程除了有限元最难可能就是风荷载了。这里面有一个计算参数叫SFL，横风向广义力功率谱，规范附录给出的计算公式与图表是不对应的，一直想不到原因，现在帖出来，大家来研究一下，到底是什么原因导致的。 有兴趣的朋友下面可以参与评论。图片有可能不清楚，可以直接下载PDF查看。

【点击下载】关于横风向广义力功率谱SFL参数的计算分析.pdf

【小编】辛展文同学 (华南理工大学)
 
分析背景：
在2012年荷载规范中，计算矩形截面横风向风荷载时，需要使用横风向广义力功率谱SFL，该值可以通过公式法计算得到，涉及建筑高度、宽度、进深和折算频率等参数，计算过程较为复杂。
同时，规范也给出H/√BD=6下的横风向广义力功率谱图表，供设计人员根据折算频率与深宽比查表确定。



2012荷载规范图

但是，在编程计算的过程中，发现根据广东规范给出的公式计算得到的横风向广义力功率谱SFL和图表查得数值不完全相符。
通过观察发现规范图表上的曲线不是很光滑，于是猜测这和绘图的精度有关。

广义力功率谱SFL计算公式：

Matlab编程分析根据规范给出的公式，使用Matlab编程绘制等值线图来验证猜想。
首先使用Matlab对一个500*500的矩阵进行运算，得到横风向广义力功率谱如下所示，可以看出它与规范中的图形比较不一样。


Matlab计算得到的图形上的曲线更圆滑，没有出现曲线断开的情况，部分曲线与规范中的变化趋势基本一致。

Matlab 500*500 矩阵数值模拟结果 

在此基础上降低精度，当进行到20*20的矩阵计算时，终于得到一个和规范相似的图形。得到图形如下：

Matlab 20*20 矩阵数值模拟结果在20*20的矩阵运算下，得到的等值线和2012荷载规范给出的基本符合，特别是在曲率变化和曲线断开的地方，都和规范图形相一致。但值得注意的是，在某些曲线上依然有细微的差别，比如A类地貌中最右边的这条曲线，规范中的是发散曲线，而Matlab绘制得到的曲线是收敛的。

这可能是由于不同软件之间的绘图处理方法不同，或由取值不均匀造成。
请大家帮忙，找找为什么会出现比较大的差异，特别是这里。 


在0.002等值线上的差异

【总结】 

（1）2012年荷载规范给出的横风向广义力功率谱SFL图表精度较低
（2）在使用查表法确定横风向广义力功率谱SFL时，在曲线不连续、斜率突变或折算频率fL1*大于0.2的位置，宜使用公式法进行校核，以得到更加准确的结果。

最后说一句：

规范也是一本书，尽信书不如无书，对于书我们一定要保有怀疑的精神，实践出真知，自勉！

==============================

这个是我的微信公众号[Dino结构笔记]：DINOSTRU
也可以通过扫二维码进行关注一下。

【Dino结构笔记】伸臂桁架的几何优化方法
结构工程分析计算中的一题多解法

最近学习【非解构】的公众号文章【超高层伸臂桁架的那些事儿】，有所心得体会，
伸臂桁架作为重要的加强层方案在结构的刚度难以满足的时会被使用。这么好的东西，
我肯定用过了，如下图所示（左侧是伸臂桁架的上部，右侧是腰桁架）。 


主要的使用场景有如下：
（1） 普通的框筒结构的刚度满足不了，相差还比较大，靠加大截面会造成不经济的情况下采用伸臂桁架，对于腰桁架效果更好
（2） 巨柱的方案，由于要把内筒的弯矩传给外筒，需要采用伸臂桁架
（3） 吊柱体系，由于要制作大空间，或外框柱间断了，局部采用伸臂桁架连接吊柱。更合理的名称也许是悬臂桁架（英文Hanger Truss）。"
（4） 无梁楼盖的框筒高层体系，因为没有内框梁与外框梁，所以需要伸臂桁架把弯矩带到外柱上，这种情况多发生于国外的高层。     


这篇笔记讨论伸臂桁架的几何优化，也就是伸臂桁架的局部布置的优化，通过以下的过程，
你会发现采用连续体的拓扑优化与基于虚功原理的桁架简化优化得出来的结果是一样了。
以下是这个方法的详细描述。

我们从伸臂桁架X支撑开始吧。

 方法一，连续体的拓扑优化 

算例：设一个连续体9mx9m（正方形比例），上端自由作用作用两个集中力P，底部两点固定支座，对该连续体进行拓扑优化，本文采用ESO方法，笔者基于OPENSEES的二次开发小程序，下载方式如下。连续体的拓扑优化的原理是简单的，经过应力分析后，把应力值最小的一部分（通过拒绝率RR进行控制）删除掉，然后进行下一步分析，到不要再删除后提高拒绝率，如此循环下去，就可以得到一个优化后的拓扑形状，该形状的特点就是在一定外力作用下，构件的每个区块的应力是饱满的，也就是把材料用尽达到优化的要求。



以下是优化的整个过程的小动画。

关于基于OPENSEES的拓扑优化的原理与程序，可以转到这个帖子去看。
[OPENSEES]实例35_基于OPENSEES拓扑优化程序开发
http://dinochen.com/article.asp?id=315  

方法二，基于虚功原理的优化 

（1） 由于伸臂桁架的轴力较大，以下实例简化了桁架构件（其实桁架的上弦下弦是压弯构件，所以应该采用梁柱单元来做的），建立ETABSS模型如下图所示。中间节点的四个杆件要划分单元，生成S2K文件后备用。计算草图如图所示。

（2） 采用ETO程序把模型导入S2K文件，如下图所示，生成OPENSEES的分析代码。

（3） 打开之前开发的OPENSEES的桁架找形程序，对中间点的位置进行找形分析，找形分析程序通过迭代后（也就是不断改变中间点的X坐标，重新提交计算）得到最优解如图表所示。

 基于OPENSEES的桁架找形程序  ETABS与ETO程序的伸臂桁架的局部建模

 基于能量的原理算出最优的位置图表所得解X=2250mm,也就是X=0.25B。
最终，两个方法是【殊途同归】的得到相同的结论，正方形比例的伸臂桁架的最优解是X=0.25B。
我再采用图解法，得到这个桁架的内力图解的图如下所示，非常用趣。 

那么换一个问题，换成V形支撑形式的伸臂桁架，虚功的最优解在哪里呢？
以下图为算例的求解过程。V形支撑的最优点在中间 
这就是局部伸臂桁架的几何优化的计算方法的简单介绍，
以上就是这篇笔记的全部内容。

==============================

这个是我的微信公众号[Dino结构笔记]：DINOSTRU
也可以通过扫二维码进行关注一下。

【Dino结构笔记】结构图解法的学习与编程
 基于桁架结构图解法的找形优化方法的开始探索

这期继续学习【结构小站】那一篇关于【题目是：建筑的“形”与“力”—形由力生，以力塑形，以形驭力】的文章里面相关的内容。这文章里面介绍了图解静力法，又可以称为桁架结构图解法等。这个方法的详细原理可以参考如下图的这本书Graphical Analysis（出版于1921年）。还有一点，这里面的图解法的插图太漂亮了，想起以前的桁架钢构桥都是用这个方法来设计的，如下图所示。

 悉尼海港大桥（Sydney Harbour Bridge） 

经典力学教科书 图解分析法

说起来也是有趣，以前的结构力学书是有这部分内容的，由于现在电算程序为主要计算手段，教材也没有介绍这部分桁架图解法，可能是因为有了电算这部分的价值不高，但是最近几年的超高层的设计发展，桁架图解法可用于大跨或高层的结构优化设计当中，这不得不重新去学习一下图解法了。


基于图解法的结构优化的工程应用   

【原理】图解静力法的基本原理
1. 作为结构，所有外力与支座反力可以构成一个环的向量和，也就是所有力向量之和是为零的。
2. 每个结点的力需要平衡，也就是一个节点多个作用力（外力支座反力与相近构件的轴力）对这个节点求和也是为零的。
3. 桁架由许多三角形组成，每个三角形三个构件的内力和不一定平衡，不平衡的力为作用于相邻近的三角形的桁架上，但是全部桁架三角形的内力传递完成后最终的合力也是为零。

综合上述的原理，就通过简单的几何原理就发展出结构图解法。


【操作】图解静力法的基本操作。对于操作而言，再多的文字与描述也是很难记录细节的，所以我通过ipad录制了以下实例的图解法操作的整个过程，有兴趣的朋友可以查看B站的视频，连接如下。 
1. 求解结构的支座反力
2. 把外力与支座反力绘制出力线部分（一般在最右端），标出由外力划分出来的线段名（大写字母）
3. 把桁架划分成许多三角形并标出号码，在桁架是面域，在图解法的图中是一个节点。
4. 根据三角形内的桁架角度关系与已标出的力线段的关系，通过直接相关找出节点位置
5. 相邻三角形可以通过相邻边的角度查找节点位置。

通过以上步骤就可以绘制出图解法的“图”，然后就可以通过测量力线的长度得到最后的构件轴力了。通过以下的草图向大家讲述求解的过程，详细可以去看看动画展示。


【Bilibili】结构图解法的详细绘图方法的视频
https://www.bilibili.com/video/BV1Lf4y1U7FS


 
【编程】图解静力法的程序可视化编程这一次学习图解法是为了更好地学习结构优化的老方法与新方法，所以需要大量的绘制图解法的图，那么：
手画是不可能的，这辈子都不可能手画的，手画又容易错，肯定要编程啦，编程的感觉好多了，画图的速度又快，准确率又高，超喜欢编程的。（老梗如果懂的请在下面回复） 


 (1) 编程的实现：编程的基本原理就是参考三角函数以及图解法的基本过程， 程序的界面如下图所示。

 (2) 验算算例。之前学习了SOM的马克写的《高层建筑设计-以结构为建筑》这本书，这本书提到了很多优化的方法，建议大家学习，同在SOM的Baker写了这一篇文章《Structural optimization using graphic statics》就是介绍SOM如何把图解法应用到结构找形与优化上面来。那么我就展示一个这个软件绘制图解法的算例吧。

算例1：  

算例2：  

 
算例3：  

 算例4：  

 
(3) 后续开发。

现在已经教会电脑程序自已通过计算绘制这个图解法的图了，接下来就要学习如何通过它进行结构优化，我开发了一个节点偏移的功能，随时改变节点就可以得到图解法的图的变化，这个就是个入门级别的优化手段了。这个就是这一篇DINO结构笔记的全部内容。

==============================

这个是我的微信公众号[Dino结构笔记]：DINOSTRU

【Dino结构笔记】桁架结构找形的优化方法
 基于OPENSEES的桁架找形程序的编写

【计算文件内容】相关的ETABS/SAP2000/OPENSEES模型的下载

最近结构圈多了好多优秀的公众号，推荐一下这几个，如【istructure】，【结构小站】这些优秀的公众号展示了工程师对结构工程的热爱与执着，配合最近的bilibili的宣传片【后浪】，让我这个老浪要觉得更加要好好学习，最近看到结构小站的一篇关于形与力的总结，题目是：建筑的“形”与“力”—形由力生，以力塑形，以形驭力，大家可以上网搜一下。 

                                   什么时侯有空也可以试下逆吊实验法

这篇文章真的写得不错，也启发了我，特别是介绍到前面大师们的找形作品，包括“逆吊实验法”去找形等，这个应用于西班牙圣家族大教堂的方法，我到西班牙第一次看到真实的模型，真的是震惊了。这个公众号的文章让我非常想做一下找形分析的尝试，那么我们开始吧。这次的前言写的有点多，一般我说做就做的。内容过于简单，请家长陪同观看。 

                                       去西班牙一定要拍一下这张照片

简单的原理：我通过一个两端铰接的梁来介绍找形的原理。首先以一个两端铰支的梁为例，大家都知道要提高刚度需要起拱，其实起拱的方法也可以通过找形来做，通过施加重力荷载，使梁发生变形，如图所示。记录每个节点的变形位移，最后把变形量乘于以一个负倍数得到反拱量，然后用于原结构的几何设计调整。原理是很简单的，这个一个很好的工程方法，确定这个倍数是一个知识点。如果梁的找形就会变成拱，板的找形就会变成拱壳，桁架的找形就变成拱壳（曲面）桁架等等。  

静力法找形状。本文以一个桁架的算例，如图所示。采用ETABS或SAP2000进行弹性建模可以得到以下的模型。施加重力后可以得到变形图如下图所示。采用SAP2000可以把变形值导出来，这个不是反拱值，这个是每个节点位置反拱的比例，至于总反拱要用优化算法算出来。桁架没有反拱的变形大约是16mm。 

打开ETO，ETABS TO OPENSEES，把ETABS建出来的模型导出E2K，导入到ETO程序中去，通过这个程序了解到21号节点是跨中节点，我们用这个节点作为控制节点。点击生成OPENSEES的分析代码，保存为co.tcl。

 在建模的过程中保存结构的节点坐标，节点在重力荷载下的变形，这个数据作用往后每步生成OPENSEES模型的准备。作用文本格式导入自编程序用于计算。 

调试OPENSEES使程序可以运行一个算例即可，证明OPENSEES可以进行一个荷载步的静力分析。那么打开编写的小程序，小程序太简单了，主要就是写一个小程序不断改变节点计算新的OPENSEES模型。每步生成带有新节点信息的OPENSEES模型提交运算，最后输入所有桁架构件的内力，用于计算指标。还有21号节点的位移值。 以下是程序不断调用OPENSEES程序的效果。

为了得到较为优化的结果，对于桁架结构，反拱后的结构形状必须使结构满足以下公式（我就简称力流公式吧，也可以理解成能量），代表最后满应力优化后，这个值小代表结构总的用钢量少，那么这个形状就是相对合理的形状了。如果TCL语言好的同学，这个公式的计算可以采用OPENSEES的TCL命令流计算得到。



目前这个算例的反拱值是从1m到5m，通过公式找到最小值，局部细化就可以找到最小值是2.8m的反拱量，最后我们把这个反拱量导入SAP2000就可以显示出这个相对优化的形状了。是不是很简单呢。大家可以试一试。

由于OPENSEES是一个没有界面的小程序，所以计算速度还是很快的。这个反拱后的桁架是不是好熟悉，对比一下，有点像。最后的优化结果桁架的变形从16mm优化到4.277，用钢量比例从4.17优化到3.82。 

最后有没有觉得桁架的找形后的形状有点眼熟，放在以前的古建筑中就非常之和谐了。




======================================================== 

这个是我的微信公众号[Dino结构笔记]：DINOSTRU

【Dino结构笔记】

真的好久没有在这个网站上写下Dino结构笔记了，今天总结一下以前关于弹性时程分析思考。
今天主要讲我对弹性时程分析中出现的楼层剪力R形突变的看法与研究。

【PPT下载地址】：弹性时程分析顶部剪力变大的原因研究.pdf 

【计算文件内容】： attachment.rar 包括ETABS模型，地震波与反应谱，最后表格数据 

（1）一般我们进行弹性时程分析，需要进行选波，选波后进行时程分析后一定会对比楼层剪力与反应谱的区别
往往在这个过程中，我们会发现顶部的楼层剪力比反应谱得出来的剪力还要大，似乎采用CQC或SRSS这种包络法
也不能包得住，那么一般文献会说这是由于高阶振型引起的顶部放大，其实，一直我有所怀疑，原因如下。

（2）为了更好的证明我的假定，我采用数值分析的方法比较直观，如下所示，我找了一个普通的超高层结构，
对于200米的超高层，一般会开始出现这种剪力的R形突变。

（3）以下是结构的空间三阶振型，要注意，这个是前几阶空间振型，不是X方向的前几阶振型，X方向的前几阶
振型分别是T1,T4,T7，如下图所示。

（4）以下是X方向前几阶周期值在反应谱的位置，表上是每个振型对X方向变形的参与系数，可见，第一阶振型为
70%，到第14振型，参与系数为1.122%,这个数值已经表明，高阶振型对宏观的位移及剪力的影响越来越小。
理论上不存在高阶振型影响楼层剪力产生R形突变。



（5）图为X方向每个主要振型的周期所在的位置，那么接下来为了证明R形突变是前几阶振型引起的，我们采用
扰动的方式造人工波进行时程分析。


（6）我们通过对反应谱进行扰动，如图所示，然后采用造人工波软件进行人工波的制作。GM0是不加扰动的，
GM1为对第一周期位置进行扰动，如此类推。



（7）图为扰动后的造波的反应谱与目标谱的对比。
（8）经过弹性时程分析，得到X方向的楼层剪力图如下图所示，可以发现只有第一阶扰动与第二阶扰动才对楼层
剪力有宏观的影响。

（a）从第三阶开始，扰动基本不影响剪力分布，基本可以排除高阶振型对剪力R形突变的影响
（b）第一二阶的扰动，对总剪力影响较大，第一阶扰动并不会引起R形突变
（c）第二振型是关键 ，第二振型开始出现R形突变
（d）为了进一步研究，我们增加对第二阶的扰动



（9）为了进一步证明第二周期对楼层剪力的R形突变的影响，我们增加对第二周期的扰动，操作如下：

（a）降低第一阶附近的地震影响系数
（b）增加对第二阶周期附件的扰动

在天然波的选择中，这的确是一个普遍现象，在周期大于4S的情况下，大部分地震波的影响系数小于规范值
为满足总剪力的需要，很多满足条件的地震波会在1~3s的区间有较大的影响系数值 。

目前的结论：
（1）高阶振型并不会影响剪力的分布，特别不会引起R形突变
（2）通过修改反应谱造人工波，证明剪力的R形突变是第二振型与反应谱的形状引起的。
（3）天然波产生的基底剪力不足，经常会第一不够第二来凑，因此产生R形突变
（4）本文只用了一个算例，提供数值如附件所示，大家可以测试与测算
（5）本文基于数值分析得到小结，如果有大神可以采用解析公式法进行证明就更好了

本人才疏学浅，以上只是学习小结，以供参考，如有不对，请指正。

如果你们对这种文章感兴趣，请留言并关注我的微信公众号 DINOSTRU

DINO结构分析笔记（1）: 时程分析选取地震波的问题

[下载实例PDF]：点击下载此文件：Dino笔记(1)时程分析选取地震波的问题

更新时间：2014年05月02日，编辑助手：林哲

简介：DINO结构分析笔记是一个专栏，每个月更新一次，主要讲述结构分析及设计当中一些简单的理论及软件的应用，主要包括以下内容：
（1）结构分析背后的理论
（2）软件程序的实现方法
（3）软件的介绍及技巧
（4）应用于实际工程的实例
（5）介绍一些自编小程序
每次推出的专题，采用通俗易懂的表达方法介绍理论与程序实现方法，不会过多地使用数字与公式，尽量让更多的工程师去明白其中的原理。
那么，如何合理地选取地震波呢？首先我们要确定地震波的来源： PEER的地震波数据库提供了大量强震记录，这些地震波涵盖了各类场地、经过基线调整、严格复核，最重要的是该数据库提供了强大的搜索功能——只要输入目标反应谱，系统会搜索到一定样本数量的地震波，用以匹配目标反应谱，使得样本的平均反应谱与目标反应谱在一定周期范围内的均方差最小。在这个基础上，再通过ETABS程序进行批量计算，选出基底剪力满足规范要求的地震波就可以了。总的来说，选择地震波的思路就是“先选后试”。这里的“地震波”指的是天然波，人工波不需要选取，由安评报告获得，如果是科研分析，可以自已拟合生成人工波。进行时程分析时，两组选波方法：（1）3组波，1组人工波+2组天然波，（2）7组波，2组人工波+5组天然波。要注意的是，采用双向地震波时，人工波的两个方向的组合系数是1.0+0.85，而天然波只对主波进行倍数A的放大，使其峰值达到规范规定小震、中震或大震的峰值，而次方向波也按倍数A放大，倍数A不再采用0.85进行调整。

最近有空时看看最优化的书《工程优化方法》（陈卫东），之前没有多大的接触，在数学上属于极值求解方法。原来除了梯度法以外还有很多种直接方法的，所谓直接方法就是不用求梯度，求梯度对于有限元来说是很难的一件事，计算量太大，倒不如用直接法算了。以下介绍一下最近学的一种单纯形法（simplex method），最直接法中最基本的方法。计算流程图如图所示(页71)。

原理简单来说就是在超空间中（多于三维的话），假定超多边形，检查棱点是否可行且最小，通过不断的放大缩小使棱点碰到最小值。接下来通过DELPHI程序编个小程序，算一算经典的一些函数就知道了。DELPHI对于我来说相当于MATLAB，经常做一些数值试验。最近老觉得自已的抽象思维越来越不行了，没有图形是想像不出数学或力学的一些东西，也许是过分依赖电脑图形，抽象思维还是得煅炼一下。作为工程研究为主的我，面向抽象的数学，当然就是采用图形化把数学具体化啦，这个工作我发现国内的数学家很少去做，我经常上网看得的都是国外的学者用图形学来把抽象的数学表述出来，如分形，拓展，纤维束等。接下来写几个算例在博客上面，这些算例都是《工程优化方法》书上摘下来，有些改了一些，他的习题是没有答案的，这里也提供一些答案吧。已经是博士生了，但是我还是经常会买些奇奇怪怪的数学的书，做做后面的习题，当然啦，这么复杂的习题当然是用电脑编程来做了。

题目：用单纯形法求无约束方法：

编制单纯形法程序求解过程如下：
计算公式：RESULT:=x[0]*x[0]+x[1]*x[1]-x[0]*x[1]-3*x[0] ;


答案： x1=2.0，x2=1.0；
采用图形学技术，将三维曲面画出来，再用梯度颜色一画，查看一下极值，检验结果（这个图形小程序就是来检验结果的，因为习题没有答案）：


图形界面，黑色为等值线，通过颜色可以明显看出形状，最值为x=2,y=1,z=-3，数学其实很美的！


这是一张网格图与带颜色的三维图，通过空间去看它的形状。
接下来搞一个工程问题：
有一个简单的支架结构，由两斜柱组成，重力荷载为2P=20kN，管截面壁厚T为10mm，高度H与管直径D为设计变量，要求强度与稳定满足，如公式所示，结构要求最轻，即质量最小。


这是一个二变量的优化问题，可以用图形来显示，由于存在边界，所以将边界外的值变大就行了，然后就进行单纯形优化。优化结果为H=1000x1.0=1000mm，D=100x1.6377=163.7mm，质量为1720kg。图形表达如图所示。

优化结果为H=1000x1.0=1000mm，D=100x1.6377=163.7mm，质量为1720kg。图形表达如图所示。

单纯形法也适用于多变量，由于不用求偏导，适用于有限元分析，复形法从单纯形法发展而来，具有更好的适用性。接下来学习更多的优化方法，以后再分享啦。


凹下去部分就是可行解，平台部分是不可行解，也就是失稳与超强的，任务就是找到可行区的最低点。图形的表达可以让数学不那么深奥，不那么抽象，数学不仅仅是符号的运算，........待续。

最近埋头于研究宏观单元与有限元程序的编程，长时间没有写博客，对不住各位网友啦！在研究生几年写了几个自已比较满意的前后处理程序，解决的很多工程问题，但是学术价值不高，只是学习的有限元软件的应用。后来，我明白的一个道理，学习有限元最好的方法之一就是编程把它实现，而且，在长期的学习OPENSEES与PERFORM-3D的过程中发现，看到国内外这个方面的差异。“国人应当自强！”，OPENSEES发展迅速，这也让人担忧，因为国内外这个差距会越来越大。于是我想我的博士论文的一个内容就是开发一个属于自已的有限元平台MESAP（Macro Element Structure Analysis Program），译意为：基于宏观单元的结构分析程序，“ME”又表达了是“自已”的东西，双重含义。
软件基于类来编写。可不断扩展单元库与材料本构库，配有图形化的前后处理，图形界面如图所示。算例结果如图所示。


现在的单元库已经有三种梁柱单元基于柔度法纤维单元，基于刚度法法纤维单元及基于柔度法的塑性铰单元。MESAP还具有两种剪力墙单元模型：三弹簧模型（TVLEM）与多竖向弹簧模型（MVLEM），MESAP的材料库具有多种材料模型，包括Kent-Park混凝土本构模型、Pinto钢筋本构模型，其滞回曲线与OpenSEES的对比如图所示，还有一些简单的塑性铰本构如修正Clough模型等。

近期对剪力墙与梁柱等宏单元进行充分的研究与对比，将写成系统的宏观单元研究论文发表至抗震相关的期刊，希望能发表成功。国内研究宏观单元现阶段大部分还是基于DRAIN，OPENSEES，CANNY，PERFORM-3D的等软件，自主研发的软件比较少，希望我在基于MESAP平台可以开发出更好更新的宏观单元，说实话，在OPENSEES里面开发单元实在太难了，在MESAP，自已的地盘很多约定简单，所以开发单元容易呀，呵呵。不过OPENSEES实在是太好的，还是有很多地方值得学习！
开发MESAP还需要发挥我开发前后处理的功夫，把界面与功能做完善，是目前我比OPENSEES先进的其中一点，呵呵（自吹自擂了），举个例子，在OPENSEES进行低周往复，位移控制流非常难写，我在MESAP就把位移控制加载所有重要因素写在一个界面里，只要这个界面输入完整，就可以完成一系列的位移加载，改动的话也很容易，这是我自已最喜欢的设计之一，界面如图所示。

继续努力，开发MESAP，这是第一阶段，未来的路还很长！！
期待MESAP的论文发表！！