PFC导入地震波加载

这串个案例的核心其实有两个，一个是如何从外部读取数据，一个是如何将加速度转化为速度给墙体。

一、外部导入地震波

    我们可以在网页上搜索经典的地震波曲线，比如下面这个就是百度搜到的（[转载]常用地震波数据汇总（72448地震波-hach.txt）_THheat_新浪博客 (sina.com.cn)）：

    其实我们仔细看一下时间的数据是一个等差数列，这是因为采样的时候我们都是每隔相等的时间采集一个加速度值。所以第一列数据我们是不需要的，我们需要第二列数据的加速度值，可以copy下来存储为一个txt文件。如下

把这个txt文件放入我们的项目目录中，后面进行调用。

然后我们看一下导入文件的代码：

    这个其实是比较简单的，首先我们需要用一个数组(array)来作为容器存储数据，这个容器的大小肯定是不能比你读取的数据要小的。

    第二行是打开文件的命令，第二个参数0的意思是这个文件是用来读取的，第三个参数1的意思是编码格式，对应的ASCII码，文本格式的文件读取都可以用这个，二进制文件用2，fish文件格式用0。

    第三行是读取文件了，这个100是读取100个数据。

    读取完后注意要关闭文件。

二、施加地震波

    加速度对时间积分可以得到速度，这个积分的数值实现研一的数值分析课都学过。我们就简单的认为某个时间区间内加速度是恒定的，那么：

V2=V1 A*detaT

这里V2为当前时间对应的速度，V1为上一个时间对应的速度，A为加速度，detaT为这两个时间点间的时间差。

    首先我们可以使用跳跃结构来得到当前时间对应的加速度，即addVel，然后再每个cycle中叠加速度，对于一个cycle来说，时间的差值就是timestep。需要注意的是array里面存放的数据是string，需要用float来转化为浮点数。

[wp=wall.find(1)][count=1][time_record=mech.age-100][addVel=0]def updataVel    whilestepping    time=mech.age    if time>time_record 0.01 then        addVel=jiasudu(count)        time_record=time        count =1    endif    ; V2=V1 detaT*a    wall.vel.y(wp)= wall.vel.y(wp) float(addVel)*global.timestepend

我们可以监测一下速度：

def jiance    whilestepping    wallVel= wall.vel.y(wp)end

history deletehistory id 1 @timehistory id 2 @wallVel

solve time 10

可以看一下结果，墙体的变化还是我们所预想的。

因为加速度比较小，而且我们的阻尼系数比较大，所以边坡的反应并不是很大。边坡的反应并不是这篇文章的重点，提出的这种加载地震波的方式是可以参考一下的，感兴趣的话积分方法可以采用数值分析中介绍的方法去优化。