t255.cpp

此程序为研究异形柱最不利荷载角的C++程序

//Tshape 200*500*500
#include "iostream.h"
#include "iomanip.h"
#include "fstream.h"
#include "math.h"
#include "stdio.h"  
#include "stdlib.h"
const int numdanc=30000;
const int numdans=12;
const int numfai=5000;
const float zhou=0.3;
const float pi=3.14159265;



const float E=0.002; 


const float ES=200000;
const float FZ=0.79*145.138;

const float F=30*FZ;//c30
const float dl=2;
const float ds=4*dl*dl;
const float d2l=dl*2;
//***********************
//约束混凝土区的本构关系

/*float strc(float e,float z,float fc)
{
float a,nn,l; 


nn=1.2*fc*(1-z*(-e-0.003));
l=0.24*fc;
if(e>=-0.003&&e<0) a=-1.2*fc*((-2*e)/(1.5*E)-(e*e)/(1.5*1.5*E*E));
else if(e>0) a=0;
else
{
if(nn>l)
a=-nn;
else
a=-l;
}
return(a);
}*/

float strc(float e,float z,float kc,float fc)
{
float a,nn,l; 


nn=kc*fc*(1-z*(-e-0.002*kc));
l=0.002*kc+0.15/z;
if(e>=-0.002*kc&&e<0) a=-kc*fc*((-2*e)/(kc*E)-(e*e)/(kc*kc*E*E));
else if(e>0) a=0;
else
{
if(e>-l)
a=-nn;
else
a=0;
}
return(a);
}

 //***********************
//非约束保护层混凝土区的本构关系 

float strc1(float e,float fc)
{
  float a;      
  if(e>=-0.002&&e<=0) a=-fc*((-2*e)/E-(e*e)/(E*E));    
  else if(e>=-0.0033&&e<-0.002) a=-fc*(1-0.8*(fabs(e)-E)/E);
  else a=0; 
   return(a);  	  
 }

//***************************
//钢筋本构关系
float strs(float e,float fy)
  {
float as;

if(e<fy/ES&&e>-fy/ES)
as=ES*e;
else if(e<=-fy/ES&&e>=-0.016)
as=-fy;
else if(e<-0.016) as=-(fy+0.01*ES*(fabs(e)-0.016));
else if(e>=fy/ES&&e<=0.016) as=fy;
else  as=fy+0.01*ES*(e-0.016);

return(as); 
  }

//************************
//主程序
void  main()
{
 float M,N,Q,P,H,J,C;
 int cdd,ndd,mcdd,mdd,ncdd,jdd,hdd,hcdd,jcdd,qcdd,qdd;
	 float ts,fy;
int i,j,k;
int kk,kkk;
int z0;//保护层单元数

float r,ri,r1,dr;
float fai,dfai,sita,dsita,fai1,sita1;
double nc,ns1,nsc,np;
double mcx,mcy;

double msx,msy,mscx,mscy;
double msx1,msy1;
double mx,my;
float eic[numdanc],eis[numdans];
//***************************
float eicm,eism;
int q,u,a,w,aa;
char quf;
int end;
float miu;
float shi,shi1[3000];
float wucha;
int tt,km;
float mfai[numfai][6],mmax[3000],mmax1[3000];
float fac[numfai],fag[numfai];


float facg;
//***************************************
//strs1数组储存钢筋屈服时各钢筋单元的应变
//strs2数组储存第一批混凝土单元屈服时各钢筋单元的应变
//strs3数组储存极限状态个钢筋单元的应变

float strs1[numdans],strs2[numdans],strs3[numdans];

float *ang1,*ang2;
float *dist1,*dist2;
float anga,angb,dista,distb;
float angle,dist;

float ix,iy;
float cossita,sinsita;
int ik;//混凝土单元数
float danc[numdanc][2],dans[numdans][2];//混凝土钢筋单元坐标

float n,m,p,as1,as2,as3,x0,y0,d0;
float l1,l2;//箍筋长度
float rou;//配筋率

float esb,sb;//配箍率和s/d
int num1;
float fc,d,B,s,alf;//

float e50u,e50h,z;


fy=300;
k=0;
fc=14.3; 
d=16.0;//纵筋直径  
s=100;//箍筋间距 
alf=-0.9*pi;
d0=8.0;//箍筋直径

M=500;//底肢
N=200;//肢宽
Q=500;//纵肢
P=Q-N;
H=(M-N)/2;
J=(M+N)/2;
C=30;  //保护层
B=N-C*2+d0*2; //核心宽度



//单元节点
 cdd=C/d2l+0.5;
 ndd=N/d2l+0.5;
 mcdd=(M-C)/d2l+0.5;
 mdd=M/d2l+0.5;
 ncdd=(N-C)/d2l+0.5;
 jdd=J/d2l+0.5;
 hdd=H/d2l+0.5;
 hcdd=(H+C)/d2l+0.5;
 jcdd=(J-C)/d2l+0.5;
 qcdd=(Q-C)/d2l+0.5;
 qdd=Q/d2l+0.5;


ix=0,iy=0;

//***********************
//形成钢筋单元坐标数组


dans[0][0]=C+d/2;dans[0][1]=C+d/2;
dans[1][0]=C+d/2;dans[1][1]=N-C-d/2;
dans[2][0]=M-C-d/2;dans[2][1]=C+d/2;
dans[3][0]=M-C-d/2;dans[3][1]=N-C-d/2;
dans[4][0]=H+C+d/2;dans[4][1]=C+d/2;
dans[5][0]=H+C+d/2;dans[5][1]=N-C-d/2;
dans[6][0]=J-C-d/2;dans[6][1]=C+d/2;
dans[7][0]=J-C-d/2;dans[7][1]=N-C-d/2;
dans[8][0]=H+C+d/2;dans[8][1]=Q-C-d/2;
dans[9][0]=J-C-d/2;dans[9][1]=Q-C-d/2;
dans[10][0]=H+C+d/2;dans[10][1]=N+P/2;
dans[11][0]=J-C-d/2;dans[11][1]=N+P/2;

for(i=0;i<numdans;i++)
cout<<dans[i][0]<<" "<<dans[i][1]<<" ";
//形成混凝土单元坐标数组
//1****************
for(i=0;i<cdd;i++)
for(j=0;j<ndd;j++)
{
danc[k][0]=d2l*i+dl;
danc[k][1]=d2l*j+dl;
ix=ix+danc[k][0];
iy=iy+danc[k][1];
k++;
}
//2****************
for(i=cdd;i<mcdd;i++)
for(j=0;j<cdd;j++)
{
danc[k][0]=d2l*i+dl;
danc[k][1]=d2l*j+dl;
ix=ix+danc[k][0];
iy=iy+danc[k][1];
k++;
}
//3****************
for(i=mcdd;i<mdd;i++)
for(j=0;j<ndd;j++)
{
danc[k][0]=d2l*i+dl;
danc[k][1]=d2l*j+dl;
ix=ix+danc[k][0];
iy=iy+danc[k][1];
k++;
}
//4****************
for(i=cdd;i<hdd;i++)
for(j=ncdd;j<ndd;j++)
{
danc[k][0]=d2l*i+dl;
danc[k][1]=d2l*j+dl;
ix=ix+danc[k][0];
iy=iy+danc[k][1];
k++;
}
//5****************
for(i=jdd;i<mcdd;i++)
for(j=ncdd;j<ndd;j++)
{
danc[k][0]=d2l*i+dl;
danc[k][1]=d2l*j+dl;
ix=ix+danc[k][0];
iy=iy+danc[k][1];
k++;
}
//6****************
for(i=hdd;i<hcdd;i++)
for(j=ncdd;j<qcdd;j++)
{
danc[k][0]=d2l*i+dl;
danc[k][1]=d2l*j+dl;
ix=ix+danc[k][0];
iy=iy+danc[k][1];
k++;
}
//7****************
for(i=jcdd;i<jdd;i++)
for(j=ncdd;j<qcdd;j++)
{
danc[k][0]=d2l*i+dl;
danc[k][1]=d2l*j+dl;
ix=ix+danc[k][0];
iy=iy+danc[k][1];
k++;
}
//8****************
for(i=hdd;i<jdd;i++)
for(j=qcdd;j<qdd;j++)
{
danc[k][0]=d2l*i+dl;
danc[k][1]=d2l*j+dl;
ix=ix+danc[k][0];
iy=iy+danc[k][1];
k++;
}
//保护层的单元坐标数组形成完毕
z0=k;
//9****************
for(i=cdd;i<mcdd;i++)
for(j=cdd;j<ncdd;j++)
{
danc[k][0]=d2l*i+dl;
danc[k][1]=d2l*j+dl;
ix=ix+danc[k][0];
iy=iy+danc[k][1];
k++;
}
//10***************
for(j=ncdd;j<qcdd;j++)
for(i=hcdd;i<jcdd;i++)
{
danc[k][0]=d2l*i+dl;
danc[k][1]=d2l*j+dl;
ix=ix+danc[k][0];
iy=iy+danc[k][1];
k++;
}
 ik=k;
for(i=0;i<ik;i++)
cout<<danc[i][0]<<" "<<danc[i][1]<<" ";
//输入原始数据


ts=ik*ds;
//n=zhou*fc*ts; //  需改
n=1312130;
x0=ix*ds/ts;y0=iy*ds/ts;


as1=d*d*pi/4.0;
as2=d0*d0*pi/4.0;
as3=12*12*pi/4.0;
l1=(N-2*C)*2+(M-2*C)*2;
l2=(Q-2*C)*2+(N-2*C)*2;
p=as2*(l1+l2)/((M+Q-N-2*C)*(N-2*C)*s);
rou=(10*as1+2*as3)/ts;

sb=s/d;
esb=44200*pow(sb,-0.412)/1000000.0;
cout<<esb<<" ";

/*e50u=(3+0.003*1.2*F)/(1.2*F-1000);
e50h=0.75*p*sqrt(B/s);
z=0.5/(e50u+e50h-0.0030);*/
float fcy,kc;
fcy=0.79*30.0;//
kc=1+p*235/fcy;
z=0.5/((3+0.29*fcy)/(145*fcy-1000)+0.75*p*sqrt(B/s)-0.002*kc);
int nnnn,nnn;
int num2,num3,num4;
int conc[numdanc],stel[numdans],stelp[numdans];
km=0;


do{
	//申请动态内存
ang1=new float[numfai];
ang2=new float[numfai];
dist1=new float[numfai];
dist2=new float[numfai];
kk=1;kkk=1;
 mmax[km]=0.0;
tt=0;

end=0;
fai=50e-7;
dfai=2e-7;
sita1=alf;

r1=25;

nnn=0;

num1=0;
do{


fai=fai+dfai;

//q=1;k=0;u=0;
dsita=0;
do{
nnnn=0;