高层建筑绝对刚性基础共同作用分析程序及算例
一、绝对刚性基础与地基的共同作用分析方法
大量工程实测结果表明，高层建筑箱形基础纵向弯曲值甚小，约在0.07‰～0.33‰范围，整体刚度很大是其主要特点之一。因此，高层建筑箱（筏）基础与地基的共同作用分析可简化为绝对刚性板与地基的共同作用来处理。
绝对刚性基础受力、沉降和倾斜时，基底仍然保持为平面，且与地基紧密接触，根据变形协调条件及弹性半空间模型集中地基反力与沉降的本构关系，推得绝对刚性基础（板）的地基反力与变形的计算公式如下：
式中：
xi，yi分别为第i网格中点的坐标；
{R}为各网格集中地基反力列向量；
A、B、C分别为基础平面坐标原点沿轴的整体倾斜及沉降量；
P、Mx、My分别为板上部荷载的合力及其对轴力矩的代数和；
fij为弹性半空间地基模型的柔度系数，按下式计算：
fij=(1-μ2)Fii/πEa (i＝j)
fij=(1-μ2) /πEr (i≠j)
二、编程语言（MATLAB6.0）简介及问题的程序化过程
1、MATLAB6.0简介
MATLAB即英文Matrix Laboratory的缩写，其意思是“矩阵实验室”，最开始是专门用于矩阵计算的软件。随着MATLAB推向市场，它逐渐具备了强大的数值计算及数据图示功能。
MATLAB可用以解决各种矩阵运算、微积分计算、线性与非线性方程（组）的求解、常微分及偏微分方程（组）求解、插值及拟合、统计与优化问题。从这个意义上讲，MATLAB是真正面向21世纪的科学计算语言。
2、问题的程序实现途径
问题的解决实际上关键是求解线性方程组，首先要形成方程组的系数矩阵[A]及荷载列向量{L}，为此在直角坐标系中将基础底板分为m=s×t块，观察方程组系数矩阵的特点，可以分为如下的几块：
左上角为地基的柔度矩阵Fm×m，右下角为一3×3的O阵；设左下角为矩阵C3×m，则右上角恰为-CTm×3，在MATLAB中不难实现这些矩阵的形成及最后将各子块组合在一起的过程。
在方程组形成之后，解方程组对于MATLAB来说更是不成问题，只需将程序运行过程中所需的数据写入DATA.M文件中，启动MATLAB命令窗口，在＞＞提示符下输入DATA↙，再输入base（程序名）↙，则MATLAB将自动执行程序，很快地给出排列好的各网格中心集中地基反力(kN)矩阵、平均地基反力(kPa)矩阵，并绘制地基反力分布状况的三维图形。
由于MATLAB程序是很多过程的批处理性质的文件，所以无需分很多步骤，程序实现的过程可简单地表示如下：
启动MATLAB→读入数据文件data.m→执行计算及绘图程序base.m→显示集中及分布地基反力矩阵→绘制地基反力分布三维图形→程序运行结束
三、算例
下面结合一个算例来说明程序的运行过程：
某工程箱形基础的平面尺寸为50.1m×9.8m，在图示坐标系中分成40个网格，采用弹性半空间地基模型，变形模量取8.36Mpa，泊松比取0.5，平均地基反力取140kPa。
利用MATLAB程序base.m求解，其数据文本文件为：
miu=0.5; E=8.360 ;L=50.1; B=9.8;s=5;t=8;ex=0.0;ey=0.0;P=68.7372;
程序运行后给出如下结果：
集中反力矩阵：
1.0e+003 *
3.3380 2.5885 2.4174 2.3564 2.3564 2.4174 2.5885 3.3380
1.3254 0.8746 0.8093 0.7870 0.7870 0.8093 0.8746 1.3254
1.8148 1.2573 1.1671 1.1359 1.1359 1.1671 1.2573 1.8148
1.3254 0.8746 0.8093 0.7870 0.7870 0.8093 0.8746 1.3254
3.3380 2.5885 2.4174 2.3564 2.3564 2.4174 2.5885 3.3380
地基倾斜及沉降：
0.0000 0.0000 0.00025
单位面积的地基反力矩阵：
271.9489 210.8835 196.9467 191.9786 191.9786 196.9467 210.8835 271.9489
107.9798 71.2515 65.9348 64.1175 64.1175 65.9348 71.2515 107.9798
147.8546 102.4330 95.0872 92.5427 92.5427 95.0872 102.4330 147.8546
107.9798 71.2515 65.9348 64.1175 64.1175 65.9348 71.2515 107.9798
271.9489 210.8835 196.9467 191.9786 191.9786 196.9467 210.8835 271.9489
根据上述结果，利用MATLAB的绘图功能，可以绘制出地基反力的分布示意图：
由图中可以看出：
角点处的反力最大，其次为边上的点，内部点的反力最小。这与所采用的地基模型——弹性半空间地基模型是相一致的，直观地说明了计算程序的正确性。
附：程序base.m及简单说明
%绝对刚性基础共同作用分析计算程序
data; %调用数据文件
m=s.*t;
x=zeros(1,m); %形成横坐标向量
for j=1:t
x(1,j)=0.5*L./t+(j-1).*L/t;
end
for j=t+1:m
x(1,j)=x(1,j-t);
end
y=zeros(1,m); %形成纵坐标向量
for j=1:m
y(1,j)=0.5*B./s;
end
for j=t+1:m
for i=1:s
if j>i.*t
y(1,j)=(2.*i+1).*0.5*B./s;
end
end
end
r=zeros(m,m); %计算i、j网格中心距
for i=1:m
for j=1:m
r(i,j)=sqrt((x(1,i)-x(1,j)).^2+(y(1,i)-y(1,j)).^2);
r(i,i)=1;
end
end
K=(1-miu.^2)./(E.*pi); %形成柔度矩阵F
a=B./s;
b=L./t;
w=a./b;
Fii=2.*w.*(log(1./w)+(1./w).*log(w+sqrt(w.^2+1))+log(1+sqrt(w.^2+1)));
F=zeros(m,m);
for i=1:m
for j=1:m
F(i,j)=K./r(i,j);
F(i,i)=K.*Fii./a;
end
end
C=ones(3,m); %形成系数矩阵左下角子块
for j=1:m
C(1,j)=x(1,j);
C(2,j)=y(1,j);
End
O=zeros(3,3); %形成系数矩阵右下角O子块
A=[F,-C';C,O]; %形成系数矩阵
My=P.*(ex+0.5.*L); %形成荷载列向量
Mx=P.*(ey+0.5.*B);
PL=zeros(m,1);
Lo1=[My,Mx,P]';
Lo=[PL;L1];
q=(A\Lo)' %求解线性方程组
q1=q./((L./t).*(B./s))
surfl(q1); %绘制地基反力分布图
