code.c

有限元分析源代码

            free((char *) d[3].q->dimen);
            free((char *) d[2].q->dimen);
            free((char *) d[1].q->dimen);
            free((char *) d[0].q->dimen);
        }
        free((char *) d[3].q);
        free((char *) d[2].q);
        free((char *) d[1].q);
        free((char *) d[0].q);
      break;

      case 5:
        d[4]  = Pop();
        d[3]  = Pop();
        d[2]  = Pop();
        d[1]  = Pop();
        d[0]  = Pop();
        d3.q = (*(QUANTITY * (*)()) *pc)(d[0].q, d[1].q,d[2].q, d[3].q, d[4].q, (QUANTITY *) NULL);
        if( UNITS_SWITCH == ON ) {
            free((char *) d[4].q->dimen->units_name);
            free((char *) d[3].q->dimen->units_name);
            free((char *) d[2].q->dimen->units_name);
            free((char *) d[1].q->dimen->units_name);
            free((char *) d[0].q->dimen->units_name);
            free((char *) d[4].q->dimen);
            free((char *) d[3].q->dimen);
            free((char *) d[2].q->dimen);
            free((char *) d[1].q->dimen);
            free((char *) d[0].q->dimen);
        }
        free((char *) d[4].q);
        free((char *) d[3].q);
        free((char *) d[2].q);
        free((char *) d[1].q);
        free((char *) d[0].q);
      break;
    }

     pc = pc+1;
     Push(d3);
}

int BltinVar_String()
{
DATUM d1, d[5], d3;
int no_args, i;

    d1 = Pop();
    no_args = (int) d1.sym;

    switch(no_args) {
      case 1:
        d[0]  = Pop();
        d3.q = (*(QUANTITY * (*)()) *pc)((char *) d[0].sym, (char *) NULL);
        free((char *) d[0].sym);
      break;
      
      case 2:
        d[1]  = Pop();
        d[0]  = Pop();
        d3.q = (*(QUANTITY * (*)()) *pc)((char *) d[0].sym, (char *) d[1].sym, (char *) NULL);
        free((char *) d[1].sym);
        free((char *) d[0].sym);
      break;

      case 3:
        d[2]  = Pop();
        d[1]  = Pop();
        d[0]  = Pop();
        d3.q = (*(QUANTITY * (*)()) *pc)((char *) d[0].sym, (char *) d[1].sym,(char *) d[2].sym, (char *) NULL);
        free((char *) d[2].sym);
        free((char *) d[1].sym);
        free((char *) d[0].sym);
      break;

      case 4:
        d[3]  = Pop();
        d[2]  = Pop();
        d[1]  = Pop();
        d[0]  = Pop();
        d3.q = (*(QUANTITY * (*)()) *pc)((char *) d[0].sym, (char *) d[1].sym, (char *) d[2].sym,
                            (char *) d[3].sym, (char *) NULL);
        free((char *) d[3].sym);
        free((char *) d[2].sym);
        free((char *) d[1].sym);
        free((char *) d[0].sym);
      break;

      case 5:
        d[4]  = Pop();
        d[3]  = Pop();
        d[2]  = Pop();
        d[1]  = Pop();
        d[0]  = Pop();
        d3.q = (*(QUANTITY * (*)()) *pc)((char *) d[0].sym,(char *) d[1].sym, (char*) d[2].sym,
                            (char *) d[3].sym,(char *) d[4].sym, (char *) NULL);
        free((char *) d[4].sym);
        free((char *) d[3].sym);
        free((char *) d[2].sym);
        free((char *) d[1].sym);
        free((char *) d[0].sym);
      break;
    }
     pc = pc+1;
     Push(d3);
}


/* 
 * ------------------------------------
 * Functions for Engineering Dimensions
 * ------------------------------------
 */ 

/* Push_Dimension() pushes a data type Datum onto the 
   interpreter stack. Increment the pointer to prog. */

int Push_Dimension()
{
DATUM    d;
int length;

    d.q        = (QUANTITY *) MyCalloc(1,sizeof(QUANTITY));
    switch(CheckUnits()) {
      case ON:
         d.q->dimen = (DIMENSIONS *) MyCalloc(1,sizeof(DIMENSIONS));
         UnitsCopy( d.q->dimen, ((SYMBOL *) *pc)->u.q->dimen );
         break;
      case OFF:
         d.q->dimen = (DIMENSIONS *)NULL; 
         break;
      default:          
         FatalError("In Push_Dimension() : CheckUnits is not ON or OFF",(char *)NULL);
         break;
    }
    pc = pc + 1;
    Push(d);
}

int Push_Dimensionless()
{
DATUM d;

    d.q        = (QUANTITY *) MyCalloc(1,sizeof(QUANTITY));
    if(CheckUnits() == ON ) {
        d.q->dimen = (DIMENSIONS *) MyCalloc(1,sizeof(DIMENSIONS));
        d.q->dimen->units_name   = (char *)NULL;
        d.q->dimen->units_type   = SI_US; /* A NUMBER is default as SI_US for now */
        d.q->dimen->scale_factor = 1.0;
        d.q->dimen->length_expnt = 0.0;
        d.q->dimen->mass_expnt   = 0.0;
        d.q->dimen->time_expnt   = 0.0;
        d.q->dimen->temp_expnt   = 0.0;
    }
    else
        d.q->dimen = (DIMENSIONS *)NULL;

    Push(d);
}

int Dimension_Eval()
{
DATUM    d1, d2;
int      length;
int   UnitsType;

    d1 = Pop();   /* quantity number, d1.q */
    d2 = Pop();   /* dimensions     , d2.q */

    switch(CheckUnits()) {
      case ON:
         UnitsType = CheckUnitsType();
         if( (d2.q->dimen->units_type != UnitsType) &&
             (d2.q->dimen->units_type != SI_US    ) &&
             (!strcmp(d2.q->dimen->units_name, "deg_C") ||
              !strcmp(d2.q->dimen->units_name, "deg_F")) ) {

            d2.q->value = d1.q->value*d2.q->dimen->scale_factor;
            d2.q->value = ConvertTempUnits(d2.q->dimen->units_name, d2.q->value, UnitsType);
            d1.q->value = d2.q->value;

         } else
            d1.q->value = d1.q->value*d2.q->dimen->scale_factor;

         UnitsSimplify( d2.q->dimen );
         UnitsCopy( d1.q->dimen, d2.q->dimen );

         free((char *) d2.q->dimen->units_name);
         free((char *) d2.q->dimen);
         free((char *) d2.q);
         Push(d1);
         break;
      case OFF:
         free((char *) d2.q);
         Push(d1);
         break;
      default:          
         FatalError("In Dimension_Eval() : CheckUnits is not ON or OFF",(char *)NULL);
         break;
    }
}

int Dimension_Mult()
{
DATUM d1, d2, d3;
int       length;

   if(CheckUnits() == OFF)
      FatalError("You should set units on to use this function",
                 "Fail to set units on in Dimension_Mult()",(char *)NULL);

   d2 = Pop();
   d1 = Pop();
   d3.q = (QUANTITY *) MyCalloc(1,sizeof(QUANTITY));
   d3.q->value = d1.q->value * d2.q->value;
   d3.q->dimen = UnitsMult( d1.q->dimen, d2.q->dimen );
   free((char *) d1.q->dimen->units_name);
   free((char *) d1.q->dimen);
   free((char *) d1.q);
   free((char *) d2.q->dimen->units_name);
   free((char *) d2.q->dimen);
   free((char *) d2.q);
   Push(d3);
}

int Dimension_Div()
{
DATUM d1, d2, d3;
int       length;

   if(CheckUnits() == OFF)
      FatalError("You should set units on to use this function",
                 "Fail to set units on Dimension_Div()",(char *)NULL);

   d2 = Pop();
   d1 = Pop();
   d3.q = (QUANTITY *) MyCalloc(1,sizeof(QUANTITY));
   d3.q->value = d1.q->value/d2.q->value;
   d3.q->dimen = UnitsDiv( d1.q->dimen, d2.q->dimen, NO );
   free((char *) d1.q->dimen->units_name);
   free((char *) d1.q->dimen);
   free((char *) d1.q);
   free((char *) d2.q->dimen->units_name);
   free((char *) d2.q->dimen);
   free((char *) d2.q);
   Push(d3);
}

int Dimension_Power()
{
DATUM d1, d2, d3;
double      temp;
int       length;

   if(CheckUnits() == OFF)
      FatalError("You should set units on to use this function",
                 "Fail to set units on in Dimension_Power()",(char *)NULL);

   d2 = Pop();  /* power value */
   d1 = Pop();  /* dimension */
   d3.q = (QUANTITY *) MyCalloc(1,sizeof(QUANTITY));
   d3.q->value = pow(d1.q->value, d2.q->value);
   d3.q->dimen = UnitsPower( d1.q->dimen, d2.q->value, NO );
   free((char *) d1.q->dimen->units_name);
   free((char *) d1.q->dimen);
   free((char *) d1.q);
   free((char *) d2.q->dimen->units_name);
   free((char *) d2.q->dimen);
   free((char *) d2.q);
   Push(d3);
}

/* 
 *  ------------------------
 *  Builtin Matrix Functions 
 *  ------------------------
 */ 

int Bltin_Matrix_Add()
{
DATUM d1, d2, d3;

   d2 = Pop();
   d1 = Pop();
   d1.m = MatrixAddReplace(d1.m, d2.m);
   MatrixFree(d2.m);
   Push(d1);
}

int Bltin_Matrix_Sub()
{
DATUM d1, d2, d3;

   d2 = Pop();
   d1 = Pop();
   d1.m = MatrixSubReplace(d1.m, d2.m);
   MatrixFree(d2.m);
   Push(d1);
}

int Bltin_Matrix_Mult()
{
DATUM d1, d2, d3;

   d2 = Pop();
   d1 = Pop();
   d3.m = MatrixMult(d1.m, d2.m);
   MatrixFree( d1.m );
   MatrixFree( d2.m );
   Push(d3);
}

int Bltin_Matrix_Quan_Mult()
{
DATUM d1, d2, d3;

   d1 = Pop();      /* quantity */
   d2 = Pop();      /* matrix   */
   d3.m = MatrixQuanMult(d1.q,d2.m);
   MatrixFree(d2.m);
   if( CheckUnits() == ON ) {
     free((char *) d1.q->dimen->units_name);
     free((char *) d1.q->dimen);
   }
   free((char *) d1.q);
   Push(d3);
}

int Bltin_Quan_Matrix_Mult() 
{
DATUM d1, d2, d3;

   d2 = Pop();      /* matrix   */
   d1 = Pop();      /* quantity */
   d3.m = MatrixQuanMult(d1.q,d2.m);
   MatrixFree(d2.m);
   if( CheckUnits() == ON ) {
     free((char *) d1.q->dimen->units_name);
     free((char *) d1.q->dimen);
   }
   free((char *) d1.q);
   Push(d3);
}

int Bltin_Matrix_Quan_Div()
{
DATUM d1, d2, d3;

   d1 = Pop();      /* quantity */
   d2 = Pop();      /* matrix   */
   d3.m = MatrixQuanDiv(d2.m, d1.q);
   MatrixFree(d2.m);
   if( CheckUnits() == ON ) {
     free((char *) d1.q->dimen->units_name);
     free((char *) d1.q->dimen);
   }
   free((char *) d1.q);
   Push(d3);
}

int Bltin_Matrix_Power()
{
DATUM d1, d2, d3;

   d1 = Pop();      /* quantity */
   d2 = Pop();      /* matrix   */
   d3.m = MatrixPower(d2.m, d1.q);
   MatrixFree(d2.m);
   if( CheckUnits() == ON ) {
     free((char *) d1.q->dimen->units_name);
     free((char *) d1.q->dimen);
   }
   free((char *) d1.q);
   Push(d3);
}

int Bltin_Matrix_Trans()
{
DATUM d1, d2, d3;

   d2   = Pop();      /* matrix  */
   d3.m = MatrixTranspose(d2.m);
   MatrixFree(d2.m);
   Push(d3);
}

int Bltin_Matrix_Negate()
{
DATUM d1;

   d1 = Pop();
   d1.m = MatrixNegateReplace(d1.m);
   Push(d1);
}

int Bltin_Matrix_Affirm()
{
   Push(Pop());
}

/* 
 * --------------------------------
 * Builtin Finite Element Functions
 * --------------------------------
*/ 

/* func(quantity, matrix) , such as AddNode, FixNode, NodeLoad */
int Bltin_Node_Quant()
{
DATUM d1, d2;

     d1 = Pop();        /* matrix            , d1.m */ 
     d2 = Pop();        /* quantity : node No, d2.q */
     (*(int (*)())*pc)(d2.q->value, d1.m);

     if(CheckUnits()==ON) {
        free((char *)d2.q->dimen->units_name);
        free((char *)d2.q->dimen);
     }
     free((char *)d2.q);
     MatrixFree(d1.m);
     pc = pc + 1;
     Push(d2);
}

int Bltin_Link_Node() 
{
DATUM d1, d2, d3;

     d2 = Pop();  /* matrix of nodes which will be linked together */
     d1 = Pop();  /* matrix of restriction condition of each dof   */
     Link_Node(d1.m, d2.m);
     MatrixFree( d1.m );
     MatrixFree( d2.m );
     pc = pc+1;
     Push(d1);
}

int Bltin_Add_Elmt()
{
DATUM d1, d2, d3;

     d1 = Pop();       /* name string of elmt attribute, d1.sym */
     d2 = Pop();       /* connection matrix of elmt: nodes connected to matrix, d2.m */
     d3 = Pop();       /* elmt No., d3.q */
     Add_Elmt(d3.q->value, d2.m, (char *) d1.sym);

     if(CheckUnits()==ON) {
        free((char *)d3.q->dimen->units_name);
        free((char *)d3.q->dimen);
     }
     free((char *)d3.q);
     MatrixFree(d2.m);
     free((char *) d1.sym);
     pc = pc+1;
     Push(d3);
}

/* 
 *  --------------------------------------------
 *  Builtin Functions for FE Solution P