📄 elmt_frame3d.c
字号:
p->nodal_loads[8].value);
printf(" Mx2 = %13.5e\t My2 = %13.5e\t Mz2 = %13.5e\n",
p->nodal_loads[9].value,
p->nodal_loads[10].value,
p->nodal_loads[11].value);
printf("\n");
printf(" Axial Force : x-direction = %13.5e \n", -p->nodal_loads[0].value);
printf(" Shear Force : y-direction = %13.5e \n", p->nodal_loads[1].value);
printf(" : z-direction = %13.5e \n", p->nodal_loads[2].value);
printf("\n");
break;
default:
break;
}
}
if(isw == LOAD_MATRIX ) {
for(j = 1; j <= p->size_of_stiff; j++)
fr[j-1][0] = 0.0;
for( i=1 ; i<=p->size_of_stiff; i++ ) {
for( j=1 ; j<=p->size_of_stiff; j++ ) {
fr[i-1][0] += rot[j-1][i-1]*p->nodal_loads[j-1].value;
}
}
for(j = 1; j <= p->size_of_stiff; j++)
p->nodal_loads[j-1].value = fr[j-1][0];
}
if( UNITS_SWITCH==ON ) {
free((char *) dp_length->units_name);
free((char *) dp_length);
free((char *) dp_force->units_name);
free((char *) dp_force);
free((char *) dp_moment->units_name);
free((char *) dp_moment);
}
MatrixFreeIndirectDouble(fr, p->size_of_stiff);
MatrixFreeIndirectDouble(dlocal, p->size_of_stiff);
MatrixFreeIndirectDouble(rot, p->size_of_stiff);
break;
default:
break;
}
return(p);
}
/* =============================================================== */
/* 3D Frame Stiffness: This function returns a 12x12 rotated */
/* stiffness matrix for 3D frame elmts. */
/* */
/* Input: */
/* ea->EA ; s->matrix pointer; eiz->EI_zaxis; eiy->EI_yaxis; */
/* gj->G.J; xl->length ; rot->rotation_matrix(12x12) */
/* nst->num_of_deg_of free per element; ndf->dof's per node; */
/* =============================================================== */
#ifdef __STDC__
MATRIX *beamst3d(ARRAY *p, MATRIX *s, double ea, double eiz, double eiy, double gj, double xl, double **rot, int nst, int no_dof)
#else
MATRIX *beamst3d(p,s, ea, eiz, eiy, gj, xl, rot, nst, no_dof)
ARRAY *p;
MATRIX *s;
double **rot;
double ea, eiz,eiy,gj, xl;
int nst, no_dof;
#endif
{
int i, j, k,ii,jj,kk;
double t;
DIMENSIONS *d1, *d2;
i = no_dof + 1;
j = no_dof + 2;
k = no_dof + 3;
ii = no_dof + 4;
jj = no_dof + 5;
kk = no_dof + 6;
t = ea/xl;
#ifdef DEBUG
printf(" in beamst3d() : \n");
printf(" : EA = %lf, xl = %lf \n", ea, xl);
printf(" : t = %lf\n", t);
#endif
s->uMatrix.daa[0][0] = t;
s->uMatrix.daa[i-1][i-1] = t;
s->uMatrix.daa[0][i-1] = -t;
s->uMatrix.daa[i-1][0] = -t;
t = 12 * eiz/xl/xl/xl ;
s->uMatrix.daa[1][1] = t;
s->uMatrix.daa[j-1][j-1] = t;
s->uMatrix.daa[1][j-1] = -t;
s->uMatrix.daa[j-1][1] = -t;
t = 12 * eiy/xl/xl/xl ;
s->uMatrix.daa[2][2] = t;
s->uMatrix.daa[k-1][k-1] = t;
s->uMatrix.daa[2][k-1] = -t;
s->uMatrix.daa[k-1][2] = -t;
t = gj / xl;
s->uMatrix.daa[3][3] = s->uMatrix.daa[ii-1][ii-1] = t;
s->uMatrix.daa[3][ii-1] = s->uMatrix.daa[ii-1][3] = -t;
t = (eiy+ eiy) / xl;
s->uMatrix.daa[4][4] = s->uMatrix.daa[jj-1][jj-1] = t+t;
s->uMatrix.daa[4][jj-1] = s->uMatrix.daa[jj-1][4] = t;
t = 6 * eiy/xl/xl;
s->uMatrix.daa[2][4] = s->uMatrix.daa[4][2] = -t;
s->uMatrix.daa[2][jj-1] = s->uMatrix.daa[jj-1][2] = -t;
s->uMatrix.daa[k-1][4] = s->uMatrix.daa[4][k-1] = t;
s->uMatrix.daa[jj-1][k-1] = s->uMatrix.daa[k-1][jj-1] = t;
t = (eiz+ eiz) / xl;
s->uMatrix.daa[5][5] = s->uMatrix.daa[kk-1][kk-1] = t+t;
s->uMatrix.daa[5][kk-1] = s->uMatrix.daa[kk-1][5] = t;
t = 6 * eiz/xl/xl;
s->uMatrix.daa[1][5] = s->uMatrix.daa[5][1] = t;
s->uMatrix.daa[1][kk-1] = s->uMatrix.daa[kk-1][1] = t;
s->uMatrix.daa[5][j-1] = s->uMatrix.daa[j-1][5] = -t;
s->uMatrix.daa[j-1][kk-1] = s->uMatrix.daa[kk-1][j-1] = -t;
/* ==================== */
/* Initial units buffer */
/* ==================== */
if( CheckUnits() == ON ) {
if(UNITS_TYPE == SI) {
d1 = DefaultUnits("Pa");
d2 = DefaultUnits("m");
}
else {
d1 = DefaultUnits("psi");
d2 = DefaultUnits("in");
}
UnitsMultRep( &(s->spColUnits[0]), d1, d2 );
UnitsCopy( &(s->spColUnits[1]), &(s->spColUnits[0]) );
UnitsCopy( &(s->spColUnits[2]), &(s->spColUnits[0]) );
UnitsMultRep( &(s->spColUnits[3]), &(s->spColUnits[0]), d2 );
UnitsCopy( &(s->spColUnits[4]), &(s->spColUnits[3]) );
UnitsCopy( &(s->spColUnits[5]), &(s->spColUnits[3]) );
ZeroUnits( &(s->spRowUnits[0]) );
ZeroUnits( &(s->spRowUnits[1]) );
ZeroUnits( &(s->spRowUnits[2]) );
UnitsCopy( &(s->spRowUnits[3]), d2 );
UnitsCopy( &(s->spRowUnits[4]), d2 );
UnitsCopy( &(s->spRowUnits[5]), d2 );
UnitsCopy( &(s->spColUnits[6]), &(s->spColUnits[0]) );
UnitsCopy( &(s->spColUnits[7]), &(s->spColUnits[1]) );
UnitsCopy( &(s->spColUnits[8]), &(s->spColUnits[2]) );
UnitsCopy( &(s->spColUnits[9]), &(s->spColUnits[3]) );
UnitsCopy( &(s->spColUnits[10]), &(s->spColUnits[4]) );
UnitsCopy( &(s->spColUnits[11]), &(s->spColUnits[5]) );
UnitsCopy( &(s->spRowUnits[6]), &(s->spRowUnits[0]) );
UnitsCopy( &(s->spRowUnits[7]), &(s->spRowUnits[1]) );
UnitsCopy( &(s->spRowUnits[8]), &(s->spRowUnits[2]) );
UnitsCopy( &(s->spRowUnits[9]), &(s->spRowUnits[3]) );
UnitsCopy( &(s->spRowUnits[10]), &(s->spRowUnits[4]) );
UnitsCopy( &(s->spRowUnits[11]), &(s->spRowUnits[5]) );
free((char *) d1->units_name);
free((char *) d1);
free((char *) d2->units_name);
free((char *) d2);
}
/* ===================== */
/* Rotate stiff matrix */
/* ===================== */
s->uMatrix.daa = (double **) rotate3d(s->uMatrix.daa, rot, no_dof);
#ifdef DEBUG
printf("flag: in beamst3d() : STIFF after rotation \n");
MatrixPrintIndirectDouble( s );
#endif
return(s);
}
/* =============================================================== */
/* function beamms3d() */
/* 3D Frame Mass Matrix: This function returns a 12x12 rotated */
/* mass matrix for 3D frame elmts. */
/* */
/* Input: */
/* mbar->Density * X_area ; */
/* s->matrix pointer; */
/* xl->length ; rot->rotation_matrix(12x12) */
/* nst->num_of_deg_of free per element; ndf->dof's per node; */
/* =============================================================== */
#ifdef __STDC__
MATRIX *beamms3d(ARRAY *p, MATRIX *s, int mtype, double mbar, double xl, double rg , double **rot, int nst, int no_dof)
#else
MATRIX *beamms3d(p, s, mtype, mbar, xl, rg , rot, nst, no_dof)
ARRAY *p;
MATRIX *s;
double **rot;
double mbar, rg , xl;
int nst, no_dof, mtype;
#endif
{
int n,m, i, j, k;
double t;
DIMENSIONS *d1,*d2,*d3;
#ifdef DEBUG
printf("INFO >> In beamms3d() : no_dof= %d \n",no_dof);
#endif
switch(mtype) {
case LUMPED:
t = mbar * xl/2;
s->uMatrix.daa[0][0] = s->uMatrix.daa[1][1] = s->uMatrix.daa[2][2] = t;
s->uMatrix.daa[3][3] = t*rg*rg;
s->uMatrix.daa[4][4] = s->uMatrix.daa[5][5] = t*xl*xl/12.0;
s->uMatrix.daa[no_dof][no_dof] = s->uMatrix.daa[no_dof+1][no_dof+1] = s->uMatrix.daa[no_dof+2][no_dof+2] = t;
s->uMatrix.daa[no_dof+3][no_dof+3] = t*rg*rg;
s->uMatrix.daa[no_dof+4][no_dof+4] = s->uMatrix.daa[no_dof+5][no_dof+5] = t*xl*xl/12.0;
break;
case CONSISTENT:
/* -------------------------------------- */
/* M terms =| Mjj Mjk | */
/* | Mkj Mkk | */
/* -------------------------------------- */
/* -------------------------------------- */
/* Mjj terms */
/* -------------------------------------- */
t = mbar*xl/420;
s->uMatrix.daa[0][0] = t* 140;
s->uMatrix.daa[1][1] = s->uMatrix.daa[2][2] = t* 156;
s->uMatrix.daa[3][3] = t* 140* rg * rg ;
s->uMatrix.daa[4][4] = s->uMatrix.daa[5][5] = t* 4 * xl * xl;
s->uMatrix.daa[4][2] = s->uMatrix.daa[2][4] = - t* 22 * xl;
s->uMatrix.daa[5][1] = s->uMatrix.daa[1][5] = t* 22 * xl;
for(m=1; m <= no_dof; m++){
for(n = 1; n <= no_dof; n++){
if(n > m) /* sym terms */
s->uMatrix.daa[m-1][n-1] = s->uMatrix.daa[n-1][m-1] ;
s->uMatrix.daa[n+no_dof-1][m+no_dof-1] = s->uMatrix.daa[n-1][m-1]; /* Mkk */
}
}
/* ------------------------------------ */
/* off diagonal terms & sign adjustment */
/* ------------------------------------ */
s->uMatrix.daa[4+no_dof][2+no_dof] = s->uMatrix.daa[2+no_dof][4+no_dof] = t* 22 * xl;
s->uMatrix.daa[5+no_dof][1+no_dof] = s->uMatrix.daa[1+no_dof][5+no_dof] = - t* 22 * xl;
/*-------------------------------------------*/
/* Mkj, Mjk terms - off main diagonal terms */
/*-------------------------------------------*/
s->uMatrix.daa[no_dof][0] = t* 70;
s->uMatrix.daa[no_dof+ 1][1] = s->uMatrix.daa[no_dof+2][2] = t* 54;
s->uMatrix.daa[no_dof+ 3][3] = t* 70 * rg * rg;
s->uMatrix.daa[no_dof+ 4][4] = s->uMatrix.daa[no_dof+ 5][5] = - t* 3 * xl * xl;
s->uMatrix.daa[no_dof+ 4][2] = s->uMatrix.daa[no_dof+ 1][5] = t* 13 * xl;
s->uMatrix.daa[no_dof+ 5][1] = s->uMatrix.daa[no_dof+ 2][4] = - t* 13 * xl;
/* Symmmetric terms */
for(m = 1; m <= no_dof; m++)
for(n = 1; n <= no_dof; n++)
s->uMatrix.daa[n-1][m+no_dof-1] = s->uMatrix.daa[n+no_dof-1][m-1];
break;
default:
FatalError("In elmt_frame3d() : beamms() : Type of Mass Matrix Undefined",(char *)NULL);
break;
}
#ifdef DEBUG
printf("elmt mass[12][12]");
MatrixPrintIndirectDouble(s);
#endif
/* ---------------------- */
/* initial unit of matrix */
/* ---------------------- */
if( CheckUnits() == ON ) {
if(UNITS_TYPE == SI) {
d1 = DefaultUnits("Pa");
d2 = DefaultUnits("m");
}
else {
d1 = DefaultUnits("psi");
d2 = DefaultUnits("in");
}
d3 = DefaultUnits("sec");
UnitsMultRep( &(s->spColUnits[0]), d1, d2 );
UnitsCopy( &(s->spColUnits[1]), &(s->spColUnits[0]) );
UnitsCopy( &(s->spColUnits[2]), &(s->spColUnits[0]) );
UnitsMultRep( &(s->spColUnits[3]), &(s->spColUnits[0]), d2 );
UnitsCopy( &(s->spColUnits[4]), &(s->spColUnits[3]) );
UnitsCopy( &(s->spColUnits[5]), &(s->spColUnits[3]) );
UnitsPowerRep( &(s->spRowUnits[0]), d3, 2.0, NO );
UnitsCopy( &(s->spRowUnits[1]), &(s->spRowUnits[0]) );
UnitsCopy( &(s->spRowUnits[2]), &(s->spRowUnits[0]) );
UnitsMultRep( &(s->spRowUnits[3]), d2, &(s->spRowUnits[0]) );
UnitsCopy( &(s->spRowUnits[4]), &(s->spRowUnits[3]) );
UnitsCopy( &(s->spRowUnits[5]), &(s->spRowUnits[3]) );
UnitsCopy( &(s->spColUnits[6]), &(s->spColUnits[0]) );
UnitsCopy( &(s->spColUnits[7]), &(s->spColUnits[1]) );
UnitsCopy( &(s->spColUnits[8]), &(s->spColUnits[2]) );
UnitsCopy( &(s->spColUnits[9]), &(s->spColUnits[3]) );
UnitsCopy( &(s->spColUnits[10]), &(s->spColUnits[4]) );
UnitsCopy( &(s->spColUnits[11]), &(s->spColUnits[5]) );
UnitsCopy( &(s->spRowUnits[6]), &(s->spRowUnits[0]) );
UnitsCopy( &(s->spRowUnits[7]), &(s->spRowUnits[1]) );
UnitsCopy( &(s->spRowUnits[8]), &(s->spRowUnits[2]) );
UnitsCopy( &(s->spRowUnits[9]), &(s->spRowUnits[3]) );
UnitsCopy( &(s->spRowUnits[10]), &(s->spRowUnits[4]) );
UnitsCopy( &(s->spRowUnits[11]), &(s->spRowUnits[5]) );
free((char *) d1->units_name);
free((char *) d1);
free((char *) d2->units_name);
free((char *) d2);
free((char *) d3->units_name);
free((char *) d3);
}
/* ---------------------- */
/* rotate to global axis */
/* ---------------------- */
s->uMatrix.daa = (double **) rotate3d(s->uMatrix.daa, rot, no_dof);
return(s);
}
/* ====================================================== */
/* function tmat(rot, type, p) */
/* Transformation matrix for 2,3d element */
/* Generates a transformation matrix */
/* Input: */
/* rot : rotation matrix (12 x 12) */
/* type: elmt type */
/* p: ARRAY * storing elmt info */
/* ====================================================== */
#ifdef __STDC__
double **tmat(double **rot, int type, ARRAY *p)
#else
double **tmat(rot, type, p)
double **rot;
int type;
ARRAY *p;
#endif
{
int i,j,n,sr,ii,jj,kk,nf, ndff;
double t[4][4];
double alpha, cx,cy,cz,csa,sna,elone, el;
alpha = p->eangle.value;
sr = p->size_of_stiff;
for(i = 0; i <= 2; i++)
for(j = 0; j <= 2; j++)
t[i][j] = 0.0;
for(i = 1; i <= sr; i++)
for(n =1; n<=sr;n++)
rot[i-1][n-1] = 0.0;
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -