excn2s.c

一个用来实现偏微分方程中网格的计算库

    for (m=0; m < strlen(elem_type); m++)
      elem_blk_parms[i].elem_type[m] = 
              toupper(elem_type[m]);
    elem_blk_parms[i].elem_type[m] = '\0';

    if (strncmp(elem_blk_parms[i].elem_type,"CIRCLE",3) == 0)
    {
      elem_blk_parms[i].elem_type_val = CIRCLE;
      /* set side set node stride */
        elem_blk_parms[i].num_nodes_per_side[0] = 1;
    }
    else if (strncmp(elem_blk_parms[i].elem_type,"SPHERE",3) == 0)
    {
      elem_blk_parms[i].elem_type_val = SPHERE;
      /* set side set node stride */
        elem_blk_parms[i].num_nodes_per_side[0] = 1;
    }
    else if (strncmp(elem_blk_parms[i].elem_type,"QUAD",3) == 0)
    {
      elem_blk_parms[i].elem_type_val = QUAD;
      /* determine side set node stride */
      if (elem_blk_parms[i].num_nodes_per_elem == 4)
        elem_blk_parms[i].num_nodes_per_side[0] = 2;
      else if (elem_blk_parms[i].num_nodes_per_elem == 5)
        elem_blk_parms[i].num_nodes_per_side[0] = 2;
      else 
        elem_blk_parms[i].num_nodes_per_side[0] = 3;
    }
    else if (strncmp(elem_blk_parms[i].elem_type,"TRIANGLE",3) == 0)
    {
      elem_blk_parms[i].elem_type_val = TRIANGLE;
      /* determine side set node stride */
      if (ndim == 2)  /* 2d TRIs */
      {
        if (elem_blk_parms[i].num_nodes_per_elem == 3)
          elem_blk_parms[i].num_nodes_per_side[0] = 2;
        else 
          elem_blk_parms[i].num_nodes_per_side[0] = 3;
      }
      else if (ndim == 3)  /* 3d TRIs */
      {
        elem_blk_parms[i].elem_type_val = TRISHELL;
        elem_blk_parms[i].num_nodes_per_side[0] =
          elem_blk_parms[i].num_nodes_per_elem;
      }
    }
    else if (strncmp(elem_blk_parms[i].elem_type,"SHELL",3) == 0)
    {
      elem_blk_parms[i].elem_type_val = SHELL;
      /* determine side set node stride */
      if (elem_blk_parms[i].num_nodes_per_elem == 2)
      {
        /* 2d SHELL; same as BEAM or TRUSS or BAR */
        elem_blk_parms[i].num_nodes_per_side[0] = 2;
        elem_blk_parms[i].elem_type_val = BEAM;
      }
      else if (elem_blk_parms[i].num_nodes_per_elem == 4)
        elem_blk_parms[i].num_nodes_per_side[0] = 4;
      else
        elem_blk_parms[i].num_nodes_per_side[0] = 8;
    }
    else if (strncmp(elem_blk_parms[i].elem_type,"HEX",3) == 0)
    {
      elem_blk_parms[i].elem_type_val = HEX;
      /* determine side set node stride */
      if (elem_blk_parms[i].num_nodes_per_elem == 8)
        elem_blk_parms[i].num_nodes_per_side[0] = 4;
      else if (elem_blk_parms[i].num_nodes_per_elem == 9)
        elem_blk_parms[i].num_nodes_per_side[0] = 4;
      else if (elem_blk_parms[i].num_nodes_per_elem == 12)  /* HEXSHELL */
        elem_blk_parms[i].num_nodes_per_side[0] = 4;
      else if (elem_blk_parms[i].num_nodes_per_elem == 27)
        elem_blk_parms[i].num_nodes_per_side[0] = 9;
      else
        elem_blk_parms[i].num_nodes_per_side[0] = 8;
    }
    else if (strncmp(elem_blk_parms[i].elem_type,"TETRA",3) == 0)
    {
      elem_blk_parms[i].elem_type_val = TETRA;
      /* determine side set node stride */
      if (elem_blk_parms[i].num_nodes_per_elem == 4)
        elem_blk_parms[i].num_nodes_per_side[0] = 3;
      else if (elem_blk_parms[i].num_nodes_per_elem == 8)
        elem_blk_parms[i].num_nodes_per_side[0] = 4;
      else
        elem_blk_parms[i].num_nodes_per_side[0] = 6;
    }
    else if (strncmp(elem_blk_parms[i].elem_type,"WEDGE",3) == 0)
    {
      elem_blk_parms[i].elem_type_val = WEDGE;
      /* determine side set node stride */
      if (elem_blk_parms[i].num_nodes_per_elem == 6)
        elem_blk_parms[i].num_nodes_per_side[0] = 4;
      else
        elem_blk_parms[i].num_nodes_per_side[0] = 8;
      sprintf(errmsg,
             "Warning: WEDGE%d is assumed to have %d nodes per face",
              elem_blk_parms[i].num_nodes_per_elem,
              elem_blk_parms[i].num_nodes_per_side[0]);
      ex_err("ex_cvt_nodes_to_sides",errmsg,EX_MSG);
    }
    else if (strncmp(elem_blk_parms[i].elem_type,"PYRAMID",3) == 0)
    {
      elem_blk_parms[i].elem_type_val = PYRAMID;
      /* determine side set node stride */
      if (elem_blk_parms[i].num_nodes_per_elem == 5)
        elem_blk_parms[i].num_nodes_per_side[0] = 4;
      else
        elem_blk_parms[i].num_nodes_per_side[0] = 8;
      sprintf(errmsg,
             "Warning: PYRAMID%d is assumed to have %d nodes per face",
              elem_blk_parms[i].num_nodes_per_elem,
              elem_blk_parms[i].num_nodes_per_side[0]);
      ex_err("ex_cvt_nodes_to_sides",errmsg,EX_MSG);
    }
    else if (strncmp(elem_blk_parms[i].elem_type,"BEAM",3) == 0)
    {
      elem_blk_parms[i].elem_type_val = BEAM;
      /* determine side set node stride */
      if (elem_blk_parms[i].num_nodes_per_elem == 2)
        elem_blk_parms[i].num_nodes_per_side[0] = 2;
      else 
        elem_blk_parms[i].num_nodes_per_side[0] = 3;
    }
    else if ( (strncmp(elem_blk_parms[i].elem_type,"TRUSS",3) == 0) ||
              (strncmp(elem_blk_parms[i].elem_type,"BAR",3) == 0)  ||
              (strncmp(elem_blk_parms[i].elem_type,"EDGE",3) == 0) )
    {
      elem_blk_parms[i].elem_type_val = TRUSS;
      /* determine side set node stride */
      if (elem_blk_parms[i].num_nodes_per_elem == 2)
        elem_blk_parms[i].num_nodes_per_side[0] = 2;
      else 
        elem_blk_parms[i].num_nodes_per_side[0] = 3;
    }
    else if (strncmp(elem_blk_parms[i].elem_type,"NULL",3) == 0)
    {
      elem_blk_parms[i].elem_type_val = NULL_ELEMENT;
      /* set side set node stride */
      elem_blk_parms[i].num_nodes_per_side[0] = 0;
    }
    else
    { /* unsupported element type; no problem if no sides specified for
         this element block */
      elem_blk_parms[i].elem_type_val = UNK;
      elem_blk_parms[i].num_nodes_per_side[0] = 0;
    }
    elem_blk_parms[i].elem_blk_id = elem_blk_ids[i];    /* save id */
    elem_ctr += elem_blk_parms[i].num_elem_in_blk;
    elem_blk_parms[i].elem_ctr = elem_ctr;      /* save elem number max */
  }


  /* Allocate space for the ss element to element block parameter index array */
  if (!(ss_parm_ndx=malloc(tot_num_ss_elem*sizeof(int))))
  {
    free(elem_blk_parms);
    free(elem_blk_ids);
    free(ss_elem_ndx);
    exerrval = EX_MEMFAIL;
    sprintf(errmsg,
"Error: failed to allocate space for side set elem parms index for file id %d"
,
            exoid);
    ex_err("ex_cvt_nodes_to_sides",errmsg,exerrval);
    return (EX_FATAL);
  }


  /* Allocate space for the ss element to node list index array */
  if (!(ss_elem_node_ndx=malloc((tot_num_ss_elem+1)*sizeof(int))))
  {
    free(ss_parm_ndx);
    free(elem_blk_parms);
    free(elem_blk_ids);
    free(ss_elem_ndx);
    exerrval = EX_MEMFAIL;
    sprintf(errmsg,
"Error: failed to allocate space for side set elem to node index for file id %d",
            exoid);
    ex_err("ex_cvt_nodes_to_sides",errmsg,exerrval);
    return (EX_FATAL);
  }

/* determine if each side set has uniform element types; this will
   be used to help determine the stride through the node list
*/

  /* Allocate space for same element type flag array*/
  if (!(same_elem_type=malloc(num_side_sets*sizeof(int))))
  {
    free(ss_elem_ndx);
    exerrval = EX_MEMFAIL;
    sprintf(errmsg,
   "Error: failed to allocate space for element type flag array for file id %d",
            exoid);
    ex_err("ex_cvt_nodes_to_sides",errmsg,exerrval);
    return (EX_FATAL);
  }

  elem_ctr = num_elem_per_set[0];
  same_elem_type[0] = TRUE;
  for (i=0,k=0;i<tot_num_ss_elem;i++)
  {
    for (j=0; j<num_elem_blks; j++)
    {
      if (side_sets_elem_list[i] <= elem_blk_parms[j].elem_ctr) break;
    }

    if (i==0) {
      el_type = elem_blk_parms[j].elem_type_val;
    } 

    /* determine which side set this element is in; assign to kth side set */
    if (i >= elem_ctr) {
      elem_ctr += num_elem_per_set[++k];
      el_type = elem_blk_parms[j].elem_type_val;
      same_elem_type[k] = TRUE;
    }

    if (el_type != elem_blk_parms[j].elem_type_val) same_elem_type[k] = FALSE;

  }

/* Build side set element to node list index and side set element 
   parameter index.
*/
  node_ctr = 0;
  elem_ctr = num_elem_per_set[0];
  for (i=0,k=0;i<tot_num_ss_elem;i++)
  {
    for (j=0; j<num_elem_blks; j++)
    {
      if (side_sets_elem_list[i] <= elem_blk_parms[j].elem_ctr)
      {
        ss_parm_ndx[i] = j;     /* assign parameter block index */
        break;
      }
    }
    ss_elem_node_ndx[i] = node_ctr;     /* assign node list index */

    /* determine which side set this element is in; assign to kth side set */
    if (i >= elem_ctr) {
       /* skip over NULL side sets */
       while (num_elem_per_set[++k] == 0);
       elem_ctr += num_elem_per_set[k];
    }

    /* determine number of nodes per side */
    if (((num_nodes_per_set[k] % num_elem_per_set[k]) == 0) &&
         (same_elem_type[k])) {  /* all side set elements are same type */
       node_ctr += num_nodes_per_set[k] /num_elem_per_set[k];
    } else {
       node_ctr += elem_blk_parms[j].num_nodes_per_side[0];
    }
  }

  ss_elem_node_ndx[i] = node_ctr;       /* assign node list index */
  free(same_elem_type);
  
  /* All setup, ready to go ... */
  
  elem_ctr=0;

  for (j=0; j < tot_num_ss_elem; j++)
  {

    if (side_sets_elem_list[ss_elem_ndx[j]] > elem_ctr)
    {
      /* release connectivity array space and get next one */
      if (elem_ctr > 0)
        free(connect);

      /* Allocate space for the connectivity array for new element block */
      if (!(connect= 
             malloc(elem_blk_parms[ss_parm_ndx[ss_elem_ndx[j]]].num_elem_in_blk*
                 elem_blk_parms[ss_parm_ndx[ss_elem_ndx[j]]].num_nodes_per_elem*
                                 sizeof(int))))
      {
        exerrval = EX_MEMFAIL;
        sprintf(errmsg,
        "Error: failed to allocate space for connectivity array for file id %d",
                exoid);
        ex_err("ex_cvt_nodes_to_sides",errmsg,exerrval);
        free(ss_elem_node_ndx);
        free(ss_parm_ndx);
        free(elem_blk_parms);
        free(elem_blk_ids);
        free(ss_elem_ndx);
        return (EX_FATAL);
      }

      /* get connectivity array */
      if (ex_get_elem_conn(
                        exoid,
                        elem_blk_parms[ss_parm_ndx[ss_elem_ndx[j]]].elem_blk_id,
                        connect) == -1)
      {
        sprintf(errmsg,
       "Error: failed to get connectivity array for elem blk %d for file id %d",
                elem_blk_parms[ss_parm_ndx[ss_elem_ndx[j]]].elem_blk_id,
                exoid);
        ex_err("ex_cvt_nodes_to_sides",errmsg,exerrval);
        free(connect);
        free(ss_elem_node_ndx);
        free(ss_parm_ndx);
        free(elem_blk_parms);
        free(elem_blk_ids);
        free(ss_elem_ndx);
        return (EX_FATAL);
      }
      elem_ctr = elem_blk_parms[ss_parm_ndx[ss_elem_ndx[j]]].elem_ctr;
    }
/*  For the first node of each side in side set, using a linear search 
      (of up to num_nodes_per_elem) of the connectivity array, 
      locate the node position in the element. The first node position
      and the second node position are used with a element type specific
      table to determine the side. */

    elem_num = side_sets_elem_list[ss_elem_ndx[j]]-1;/* element number 0-based*/
    /* calculate the relative element number position in it's block*/
    elem_num_pos = elem_num - 
                  (elem_blk_parms[ss_parm_ndx[ss_elem_ndx[j]]].elem_ctr -
                   elem_blk_parms[ss_parm_ndx[ss_elem_ndx[j]]].num_elem_in_blk);
    /* calculate the beginning of the node list for this element by
         using the ss_elem_node_ndx index into the side_sets_node_index