计算力学——课程设计
withini=1; within1=n; } if(dof[m][n][0]==j) { withinj=1; within2=n; } if(withini*withinj==1) { Kframe[i][j]+=gk[m][within1][within2]; withini=0; withinj=0; } } } }
⑤void getload() //编出荷载列阵获取函数； 在成员变量中，我们定义了double gl[Nnode][3]为节点荷载列阵，在初始数据中我
们输入了该变量，故可直接获得节点的荷载列阵。
实际情况中，我们往往遇到的荷载并非直接施加于节点上，出现概率较大的往往
是梁上集中荷载、集中弯矩、均布荷载或线性分布荷载。但我们可以利用荷载等效原则，将各种荷载状况转换为节点荷载。 具体函数可见：equiload.cpp文件。
⑥ void caculate() //利用先处理法对总体刚度矩阵进
行处理，再求得各节点位移；
在求解平衡方程之前，我们先对刚度矩阵前处理。由于框架与剪力墙部分边界结
点是固定支承，它们的所有位移分量已知为零。故位移分量已知为零的结点，其结点编码均可设为“0”,且编码为“0”的结点一律不设置地址，只有位移分量未知的结点才可进入提供的内存空间。
按此原则，我们可以形成计算矩阵Ab[i][j]，然后根据高斯消去法求解方程，求
得各节点的位移。最后将各节点的位移反带到各单元的单元刚度矩阵中，便可求得节点力。
⑦ void print() //打印出计算结果力和位移
在此函数中，我们将各节点的位移与荷载输入到output.txt文件中，以便查询。
第 6 页 共 31 页
计算力学——课程设计
开始 输入原始数据 提取数据生成梁单元 读取节点自由度的的整体编码和局部编码 计算坐标转换矩阵 组装总刚矩阵 计算刚架单元单刚 组装总荷载向量zonghezaixiangli计算单元节点荷载列阵 解方程求解节点位移（高斯消去法） 计算节点力 输出相关计算结果 图3 框架模块流程图
（2） 剪力墙模块的具体流程图如图4所示。 Wall类由以下几部分组成： 1、成员变量：
const int Nmeshx=2; //定义剪力墙X方向上的划分数 const int Nmeshy=4; //定义剪力墙X方向上的划分数 double point[4][2]; //定义剪力墙结构的四个角点 double E,v; //定义弹性模量和泊松比 double gp[Nmeshy+1][Nmeshx+1][2]; //得到内部节点的坐标数组 int dof[Nmeshx*Nmeshy][8][2]; //获取自由度的整体编码和局部编码 double gk[Nmeshy][Nmeshx][8][8]; //得到单元刚度的数组 double Kwall[(Nmeshx+1)*(Nmeshy+1)*2][(Nmeshx+1)*(Nmeshy+1)*2]; //得到总体刚度的数组 double gl[(Nmeshx+1)*(Nmeshy+1)][2]; //获得节点荷载列阵 double loadwall[(Nmeshx+1)*(Nmeshy+1)*2]; //获得剪力墙的荷载列阵(单位为kN,m) double displacement[(Nmeshx+1)*(Nmeshy+1)*2]; //位移列阵(单位为m) double Fpoint[(Nmeshx+1)*(Nmeshy+1)*2][2]; //算出节点上的力(单位为kN,m)
第 7 页 共 31 页
计算力学——课程设计
定义剪力墙的角点坐标 划分网格，获得内部节点的坐标数组 计算节点自由度的的整体编码和局部编码 等参元 集成剪力墙总刚 计算四节点平面单元刚度矩阵 组装总荷载向量zonghezaixiangli解方程求解节点位移 计算单元节点荷载列阵 计算节点力 输出相关计算结果 图4 剪力墙模块流程图
3、 成员函数
①void inputdata() //用于输入相关数据，包括剪力墙四个角点的坐标、单元属性与节点荷载；
②void getpoint() //获取节点各自由度的局部编码与整体编码；
在程序中，我们定义了常变量Nmeshx与Nmeshy，即两个方向的划分数。根据网格划分数与单元节点数，我们可以确定各节点的坐标。其中，gp[i][j][0]存储横坐标，而gp[i][j][1]存储纵坐标。
③void getKelement() //获取单元刚度矩阵函数；
根据网格划分情况与单元节点数，我们可以得到各节点自由度的整体编码与局部编码。其中，dof[i][j][0]存储整体编码，而dof[i][j][1]存储局部编码。根据计算力学中介绍 平面四节点单元特性，我们可以计算得出其单元刚度矩阵，是一个8?8的矩阵。 ④void getKwall() //单元刚度合成整体刚度矩阵函数；
在输入数据中，我们定义了节点编码与自由度。其中自由度是按一定顺序输入的，故其本质就是每个节点两个自由度的定位向量。在定位向量的基础上，我们可以进行总体刚度矩阵的合成。具体思路是：每个单元刚度矩阵中的元素的行标与列标均有具体的定位向量与其对应。在合成总刚的时候，总体刚度矩阵中的行标与列阵就是定位向量，然后我们利用for循环在各个单元刚度矩阵中寻找与某定位向量对应的元素，然后相加得到最终的总刚。
关键函数如下： int withini,withinj,within1,within2; //判断单元是否包含两个作用的自由度 for(i=0;i<(Nmeshx+1)*(Nmeshy+1)*2;i++) {
第 8 页 共 31 页
计算力学——课程设计
for(int j=0;j<(Nmeshx+1)*(Nmeshy+1)*2;j++) { for(int m1=0;m1 ⑤void getload() //编出荷载列阵获取函数； 在成员变量中，我们定义了double gl[Nnode][3]为节点荷载列阵，在初始数据中我 们输入了该变量，故可直接获得节点的荷载列阵。 实际情况中，我们往往遇到的荷载并非直接施加于节点上，出现概率较大的往往 是梁上集中荷载、集中弯矩、均布荷载或线性分布荷载。但我们可以利用荷载等效原则，将各种荷载状况转换为节点荷载。 具体函数可见：equiload.cpp文件。 ⑥ void caculate() //利用先处理法对总体刚度矩阵进行处 理，再求得各节点位移； 在求解平衡方程之前，我们先对刚度矩阵前处理。由于框架与剪力墙部分边界结 点是固定支承，它们的所有位移分量已知为零。故位移分量已知为零的结点，其结点编码均可设为“0”,且编码为“0”的结点一律不设置地址，只有位移分量未知的结点才可进入提供的内存空间。 按此原则，我们可以形成计算矩阵Ab[i][j]，然后根据高斯消去法求解方程，求 得各节点的位移。最后将各节点的位移反带到各单元的单元刚度矩阵中，便可求得节点力。 ⑦ void print() //打印出计算结果力和位移 在此函数中，我们将各节点的位移与荷载输入到output.txt文件中，以便查询。 第 9 页 共 31 页 计算力学——课程设计 （3） 耦合模块具体流程 fw类由以下几部分组成： 1、 成员变量： const int Njoint=4; //定义节点数； int joint[Njoint][2]; //第一个元素表示剪力墙的整体坐标，第二个元素表示框架的整体编码 double Kframewall[(Nmeshx+1)*(Nmeshy+1)*3+Nnode*3][(Nmeshx+1)*(Nmeshy+1)*3+ Nnode*3]; //所有点的刚度矩阵 double Kwallch[(Nmeshx+1)*(Nmeshy+1)*3][(Nmeshx+1)*(Nmeshy+1)*3]; //剪力墙的变换后的刚度矩阵 double loadfw[(Nmeshx+1)*(Nmeshy+1)+Nnode][3]; //所有点的荷载 double displacementfw[(Nmeshx+1)*(Nmeshy+1)*3+Nnode*3]; //所有点的位移 2、成员函数 ①void inputdata() //读取框架与剪力墙耦合点的剪力墙节点编号与框架的节点编号； ② Kwallchange() //剪力墙总刚的扩充； 由于剪力墙模块中，我们使用的是平面四节点单元，每个节点只有两个自由度。而框架梁单元中每个节点有三个自由度。故如果对两个单元进行耦合，我们认为需对剪力墙的单元刚度矩阵进行扩充，即：每个节点由原来的两个自由度扩充为三个自由度，单元刚度矩阵也随之由原来的8?8矩阵扩充为12?12矩阵。扩充的那两个自由度对应刚度矩阵中的元素均置0。 ③getKfw() //所有点的总刚处理，对所有点进行新的编码，剪力墙原编码统一加上一个数16形成新的编码； 在这里我们重新定义了耦合后的刚度矩阵double Kframewall[(Nmeshx+1)*(Nmeshy +1)*3+Nnode*3][(Nmeshx+1)*(Nmeshy+1)*3+Nnode*3]。其数组大小是根据框架的总结点数与剪力墙的划分情况而确定的。 在框架总体刚度矩阵与剪力墙刚度矩阵耦合时，首先根据节点自由度判断是否为耦合点自由度，若是则剪力墙总体刚度矩阵与框架总体刚度矩阵中的相关项进行叠加。 for(i=0;i<(Nmeshx+1)*(Nmeshy+1)*3;i++) //判断耦合节点，从而将 剪力墙总刚与框架总刚进行集成 { for(int j=0;j<(Nmeshx+1)*(Nmeshy+1)*3;j++) { int mi1,mi2,mi3,mj1,mj2,mj3; mi1=i/3; mi2=i%3; mi3=0; mj1=j/3; mj2=j%3; mj3=0; for(int n=0;n 第 10 页 共 31 页