fe_print.c

有限元分析源代码

                              ZeroUnits(&(array->displ->spColUnits[i-1]));
                           }
                         }
                         break;
                    case 4:
                    case 8:
                        ii = eap->map_ldof_to_gdof[k-1];
                        jj = frame->node[node_no - 1].bound_id[ii-1];
                        if(jj > 0) {
                           array->displ->uMatrix.daa[k-1][i-1] = m1->uMatrix.daa[jj-1][0];
                           if( UNITS_SWITCH == ON ) {
                              UnitsCopy(&(array->displ->spRowUnits[k-1]), &(m1->spRowUnits[jj-1]));
                              ZeroUnits( &(array->displ->spColUnits[i-1]) );
                           }
                        }
                        else {
                           array->displ->uMatrix.daa[k-1][i-1]
                           = frame->node[node_no -1].disp[ii-1].value;
                           if( UNITS_SWITCH == ON ) {
                              UnitsCopy( &(array->displ->spRowUnits[k-1]),
                                frame->node[node_no -1].disp[ii-1].dimen);
                              ZeroUnits( &(array->displ->spColUnits[i-1]) );
                           }
                        }
                        k = k + 1;
                        break;
                      default:
                        break;
                   }
               }
           }
	   PRINT_STRESS = ON;
           array = elmlib(array, STRESS);
	   PRINT_STRESS = OFF;

           for( i=1 ; i <= frame->no_dof ; i++ )
              for( j=1 ; j <= frame->no_nodes_per_elmt ; j++ )
                  frame->element[elmt_no-1].rp->Forces->uMatrix.daa[i-1][j-1]
                  = array->nodal_loads[frame->no_dof*(j-1)+i-1].value;
           if( UNITS_SWITCH == ON ) {
              for( i=1 ; i <= frame->no_dof ; i++ )
                  UnitsCopy( &(frame->element[elmt_no-1].rp->Forces->spRowUnits[i-1]) , array->nodal_loads[i-1].dimen );
              for( j=1 ; j <= frame->no_nodes_per_elmt ; j++ )
                  ZeroUnits( &(frame->element[elmt_no-1].rp->Forces->spColUnits[j-1]) );
           }
       }
       
      if(m1 == (MATRIX *) NULL) { /* stress response from response structure */


         slp = lookup("InPlaneIntegPts");  /* no of integration pts in plane/surface */
         if(slp == NULL)
            iInPlaneIntegPts = 2*2;        /* 2x2 as default */
         else
            iInPlaneIntegPts = (int) slp->u.q->value;

         slp = lookup("ThicknessIntegPts"); /* no of integration pts in thickness direction */
         if(slp == NULL)
            iThicknessIntegPts = 2;        /* 2 as default */
         else
            iThicknessIntegPts = (int) slp->u.q->value;
  
         iNO_INTEG_pts = iInPlaneIntegPts*iThicknessIntegPts;

         lint = (int) 0;
         l = (int) iInPlaneIntegPts;
         for(i = 1; i<= l*l; i++) {
             xi[i-1]  = 0.0;
             eta[i-1] = 0.0;
             wg[i-1]  = 0.0;
         }

         if(l*l != lint)
            pgauss(l,&lint,xi,eta,wg);

         gamma = dVectorAlloc(iThicknessIntegPts+1);
         wt    = dVectorAlloc(iThicknessIntegPts+1);
         gauss(gamma,wt,iThicknessIntegPts);

         for(elmt_no = 1; elmt_no <= frame->no_elements; elmt_no++) {
             ep           = &frame->element[elmt_no-1];
             elmt_attr_no = ep->elmt_attr_no;
             eap          = &frame->eattr[elmt_attr_no-1];
             if(elmt_no == 1) 
                printf(" Element : %s \n Material : %s \n\n", eap->elmt_type, eap->material);

             printf("\n STRESS in  Element No  %d \n",elmt_no);
             printf(" =============================================================================================================== \n");
             printf(" Gaussion    xi       eta        gamma    Stre-xx         Stre-yy         Stre-xy         Stre-yz        Stre-zx \n");
             if(UNITS_SWITCH == OFF)
                printf("  Points \n");
             if(UNITS_SWITCH == ON)   {
                printf("  Points                                    %s             %s             %s              %s             %s \n",
                          ep->rp->stress->spRowUnits[0].units_name,
                          ep->rp->stress->spRowUnits[1].units_name,
                          ep->rp->stress->spRowUnits[2].units_name,
                          ep->rp->stress->spRowUnits[3].units_name,
                          ep->rp->stress->spRowUnits[4].units_name);
             printf(" ---------------------------------------------------------------------------------------------------------------\n \n");

                for(i = 1; i <= iThicknessIntegPts; i++)
                    for(j = 1; j <= lint; j++) {
                        k = lint*(i-1)+j;
                        printf("   %d  %10.4f %10.4f %10.4f", k, xi[j-1], eta[j-1], gamma[i]);
                        printf("\t%11.4e\t%11.4e\t%11.4e\t%11.4e\t%11.4e\n", 
                             ep->rp->stress->uMatrix.daa[0][k-1]/ep->rp->stress->spRowUnits[0].scale_factor,
                             ep->rp->stress->uMatrix.daa[1][k-1]/ep->rp->stress->spRowUnits[1].scale_factor,
                             ep->rp->stress->uMatrix.daa[2][k-1]/ep->rp->stress->spRowUnits[2].scale_factor,
                             ep->rp->stress->uMatrix.daa[3][k-1]/ep->rp->stress->spRowUnits[3].scale_factor,
                             ep->rp->stress->uMatrix.daa[4][k-1]/ep->rp->stress->spRowUnits[4].scale_factor);
                    }
             }
             if(UNITS_SWITCH == OFF)
                for(i = 1; i <= iThicknessIntegPts; i++)
                    for(j = 1; j <= lint; j++) {
                        k = lint*(i-1)+j;
                        printf("   %d  %10.4f %10.4f %10.4f", k, xi[j-1], eta[j-1], gamma[i]);
                        printf("\t%11.4e\t%11.4e\t%11.4e\t%11.4e\t%11.4e\n",
                             ep->rp->stress->uMatrix.daa[0][k-1],
                             ep->rp->stress->uMatrix.daa[1][k-1],
                             ep->rp->stress->uMatrix.daa[2][k-1],
                             ep->rp->stress->uMatrix.daa[3][k-1],
                             ep->rp->stress->uMatrix.daa[4][k-1]);
                    }
         }
      }

#ifdef DEBUG
       printf("*** Leave Print_Stress()\n");
#endif

}

#ifdef __STDC__
void Print_Displ(MATRIX *m1)
#else
void Print_Displ(m1)
MATRIX *m1;
#endif
{
int node_no, j, k, ii, jj;
int          UNITS_SWITCH;
double                 da;

#ifdef DEBUG
       printf("\n*** Enter Print_Displ()\n");
#endif

    UNITS_SWITCH = CheckUnits();

    if(frame->no_dof == 6) {
       printf("\n------------------------------------------------------------------------------------------------------------------\n");
       printf(" Node                                             Displacement\n");
       printf("  No             displ-x           displ-y           displ-z             rot-x             rot-y             rot-z\n");
       printf("------------------------------------------------------------------------------------------------------------------\n");

       if( UNITS_SWITCH == ON ) {
          printf(" units");
          for( ii=1 ; ii<=6 ; ii++ ) {
             if( ii <= 3 ) {
                if( (UNITS_TYPE==SI) || (UNITS_TYPE==SI_US) )
                   printf("                m ");
                else 
                   printf("               in ");
             }
             else
                   printf("              rad ");
          }
       }
       printf("\n");
    }

    if(frame->no_dof == 5) {
       printf("\n------------------------------------------------------------------------------------------------\n");
       printf(" Node                                   Displacement\n");
       printf("  No             displ-x           displ-y           displ-z             rot-x             rot-y\n");
       printf("------------------------------------------------------------------------------------------------\n");
       if( UNITS_SWITCH == ON ) {
          printf(" units");
          for( ii=1 ; ii<=5 ; ii++ ) {
             if( ii <= 3 ) {
                if( (UNITS_TYPE==SI) || (UNITS_TYPE==SI_US) )
                   printf("                m ");
                else 
                   printf("               in ");
             }
             else
                   printf("              rad ");
          }
       }
       printf("\n");
    }

    if(frame->no_dof == 3) {
       printf("\n------------------------------------------------------------\n");
       printf(" Node                           Displacement\n");
       printf("  No             displ-x           displ-y             rot-z\n");
       printf("------------------------------------------------------------\n");

       if( UNITS_SWITCH == ON ) {
          printf(" units");
          for( ii=1 ; ii<=3 ; ii++ ) {
             if( ii <= 2 ) {
                if( (UNITS_TYPE==SI) || (UNITS_TYPE==SI_US) )
                   printf("                m ");
                else 
                   printf("               in ");
             }
             else
                   printf("              rad ");
          }
       }
       printf("\n");
    }

    if(frame->no_dof == 2) {

       printf("\n------------------------------------------\n");
       printf(" Node                  Displacement\n");
       printf("  No             displ-x           displ-y\n");
       printf("------------------------------------------\n");
       if( UNITS_SWITCH == ON ) {
          printf(" units");
          for( ii=1 ; ii<=2 ; ii++ ) {
                if( (UNITS_TYPE==SI) || (UNITS_TYPE==SI_US) )
                   printf("                m ");
                else 
                   printf("               in ");
          }
       }
       printf("\n");
    }

 /*
  *  ===================
  *  Print displacements 
  *  ===================
  */
   
  switch( UNITS_SWITCH ) {
      case ON:
           for(node_no = 1; node_no <= frame->no_nodes; node_no++) {
           printf("%4d  ", node_no);
           for(j  = 1; j <= frame->no_dof; j++){
               jj = frame->node[node_no - 1].bound_id[j-1];
               if(jj > 0) {
                  if(m1->spRowUnits[jj-1].units_name != NULL ) {
                      UnitsTypeConvert(&(m1->spRowUnits[jj-1]), UNITS_TYPE);
                  }
                  RadUnitsSimplify( &(m1->spRowUnits[jj-1]) );
                  frame->node[node_no-1].disp[j-1].value = m1->uMatrix.daa[jj-1][0];
                  UnitsCopy( frame->node[node_no-1].disp[j-1].dimen, &(m1->spRowUnits[jj-1]) );
                  da =  MatrixContentScaleIndirectDouble(m1,jj, 1);
                  printf("     %13.5e", da);
               } else {
                  switch(UNITS_TYPE) {
                    case SI:
                     if( (frame->node[node_no-1].disp[j-1].dimen->units_name != NULL) &&
                        !strcmp( frame->node[node_no-1].disp[j-1].dimen->units_name, "deg_F"))

                        frame->node[node_no -1].disp[j-1].value
                        = ConvertTempUnits(frame->node[node_no-1].disp[j-1].dimen->units_name,
                         frame->node[node_no -1].disp[j-1].value, US);
                         break;
                    case US:
                         if((frame->node[node_no-1].disp[j-1].dimen->units_name != NULL) &&
                             !strcmp(frame->node[node_no-1].disp[j-1].dimen->units_name, "deg_C"))

                             frame->node[node_no -1].disp[j-1].value
                             = ConvertTempUnits(frame->node[node_no-1].disp[j-1].dimen->units_name,
                                                frame->node[node_no -1].disp[j-1].value, SI);
                          break;
                  }
                  RadUnitsSimplify( frame->node[node_no-1].disp[j-1].dimen );
                  printf("     %13.5e", frame->node[node_no -1].disp[j-1].value/
                        frame->node[node_no-1].disp[j-1].dimen->scale_factor);
                }
           }
           printf("\n");
         }
         break;
      case OFF:
         for(node_no = 1; node_no <= frame->no_nodes; node_no++) {
             printf("%4d  ", node_no);
             for(j  = 1; j <= frame->no_dof; j++){
                jj = frame->node[node_no - 1].bound_id[j-1];
                if(jj > 0) {
                   frame->node[node_no-1].disp[j-1].value = m1->uMatrix.daa[jj-1][0];
                   printf("     %13.5e", m1->uMatrix.daa[jj-1][0]);
                }
                else
                   printf("     %13.5e", frame->node[node_no -1].disp[j-1].value);
             }
             printf("\n");
         }
         break;
      default:
         break;
    }

#ifdef DEBUG
       printf("*** Leave Print_Displ()\n\n");
#endif

}