output1.cpp

Finite element program for mechanical problem. It can solve various problem in solid problem

#include "globprep.h"
#include "output.h"
#include <string.h>

long output(FILE *out, descrip *d)
{
  // Output data in mechtop format
  long err;
  long npropn, nprope;
  npropn = Npropn;
  nprope = Nprope;

  err = wr_nodes(out, npropn); // writing nodes, crsecs in the nodes and local systems
  fflush(out);
  if (err)
    return(7);

  err = wr_bocon(out, npropn); // writing boundary conditions
  fflush(out);
  if (err)
    return(7);

  err = wr_elements(out, nprope, d); // writing elements and code numbers on elements
  fflush(out);
  if (err)
    return(7);

  if (!(d->seq))
  {
    err = wr_globnodnum(out);
    if (err)
      return(7);
  }

//  err = wr_intpoints(out, nprope); // writing integration points of elements
                                   // with indeces of material parameters
  fflush(out);
  if (err)
    return(7);

  err = wr_materials(out); // writing material parameters
  fflush(out);
  if (err)
    return(7);

  if (Mp->stresscomp)
  {
    err = wr_auxpoint(out, nprope); // writing auxiliary points
    fflush(out);
    if (err)
      return(7);
    err = wr_ellcsys(out, nprope); // writing local coordinate systems in the auxiliary points
    fprintf(out, "\n");
    fflush(out);
    if (err)
      return(7);
  }
  if (Mp->straincomp)
  {
    err = wr_auxpoint(out, nprope); // writing auxiliary points
    fflush(out);
    if (err)
      return(7);
    err = wr_ellcsys(out, nprope); // writing local coordinate systems in the auxiliary points
    fprintf(out, "\n");
    fflush(out);
    if (err)
      return(7);
  }

  err = wr_crsecs(out);    // writing cross section parameters
  fflush(out);
  if (err)
    return(7);

  err = wr_load(out, npropn, nprope); // writing loads
  fflush(out);
  if (err)
    return(7);

  fclose(out);

  // deallocating used memory
  long i;
  for (i = 0; i < npropn; i++)
  {
    delete Boc[i];
    delete Boccn[i];
    delete Lno[i];
    delete Dlno[i];
    delete [] Nlsv[i];
    delete [] Nccn[i];
    delete Inic[i];
  }
  for (i = 0; i < nprope; i++)
  {
    delete Lel[i];
//    delete Dlel[i];
    delete [] Mel[i];
    delete [] Melid[i];
    delete [] Elsv[i];
    delete Edgl[i];
    delete Eps[i];
  }
  delete [] Nnp;
  delete [] Nep;
  delete [] Ndofn;
  delete [] Boc;
  delete [] Boccn;
  delete [] Lno;
  delete [] Dlno;
  delete [] Lel;
//  delete [] Dlel;
  delete [] Tel;
  delete [] Sel;
  delete [] Mel;
  delete [] Melid;
  delete [] Nmel;
  delete [] Ncs;
  delete [] Ncsi;
  delete [] Ecs;
  delete [] Ecsi;
  delete [] Nlsv;
  delete [] Elsv;
  delete [] Nccn;
  delete [] Edgl;
  delete Dbmat;
  delete Dbcrs;
  delete [] Tndofn;
  delete [] Tet;
  delete [] Nmu;
  delete [] Ncu;
  delete [] Inic;
  delete [] Enip;
  delete [] Edgl;
  delete [] Eep;
  delete [] Eps;
  delete Top;
  return (0);
}

long wr_nodes(FILE *out, long nprop)
/**<
 Function writes section with description of nodes to the text file given by out.

 @param out - pointer to the opened text file, where the data will be written
 @param nprop - number of entries in the input property file

 @retval 0 - if succes
 @retval 1 - if node hasn't assigned ndof
 @retval 2 - if node hasn't assigned cross-section
*/
{
  long i, j, k, l, *ic, ndf, ncf, nlsf, rindex;

  Tndofn =  new long[Top->nodes];
  memset(Tndofn, 0, sizeof(*Tndofn)*Top->nodes);
  ic = new long[Top->nodes];

  fprintf(out, "%ld\n", Top->nodes);
  for (i = 0; i < Top->nodes; i++)
  {
    for (j = 0; j < nprop; j++)
    // loop over all properties from input
    {
      if (Top->N[i].property == Nnp[j])
      // test if number of property matchs number from input
      {
        if (Ncs[j])
        // test if any cross-section is assigned
        {
          if ((ic[i] = Dbcrs->search_crs(Ncs[j], Ncsi[j]-1)) < 0) // searching cs in the cs database
          {
            fprintf(stderr,
                "\nError - cross-section with number %d and index %ld in node %ld not found\n", Ncs[j], Ncsi[j], i+1);
            delete [] ic;
            return(2);
          }
          else
          {
            if (Dbcrs->crs[ic[i]].propu[Ncsi[j]-1] == 0)
              Dbcrs->crs[ic[i]].npropu++; // indicator whether this cs type is used
            if (Ncu[ic[i]] == 0)
              Ncu[ic[i]]++;   // indicator whether cs type from cs database is used in the given task
            Dbcrs->crs[ic[i]].propu[Ncsi[j]-1]++; // indicator whether this index of cs is used
          }
          break;
        }
      }
    }
/*    if (ncf == 0)
    {
      fprintf(stderr, "\nError - node number %ld hasn't assigned cross-section\n", i+1);
      delete [] ic;
      return (2);
    }*/
  }
  Csrin = 0;
  for (j = 0; j < Dbcrs->numt; j++)
  {
    if (Dbcrs->crs[j].npropu)
    {
      rindex = 1;
      for (k = 0; k < Dbcrs->crs[j].nprop; k++)
      {
        if (Dbcrs->crs[j].propu[k])
        {
          Dbcrs->crs[j].ridx[k] = rindex;
          rindex++;
        }
      }
    }
  }
  Csrin = rindex;
  for (i = 0; i < Top->nodes; i++)
  {
    // writing node number and coordinates
    fprintf(out, "%6ld % e % e % e", i+1, Top->N[i].x, Top->N[i].y, Top->N[i].z);
    for (j = 0, ndf = 0; j < nprop; j++)
    // loop over all properties from input
    {
      if (Top->N[i].property == Nnp[j])
      // test if number of property matchs number from input
      {
        if (Ndofn[j])
        // test if node has assigned ndofs
        {
          fprintf(out, " %ld", Ndofn[j]);
          Tndofn[i] = Ndofn[j];
          ndf = 1;
          break;
        }
      }
    }
    if (ndf == 0)
    {
      fprintf(stderr, "\nError - node number %ld hasn't assigned ndofs\n", i+1);
      delete [] ic;
      return (1);
    }
    for (j = 0, ncf = 0; j < nprop; j++)
    // loop over all properties from input
    {
      if (Top->N[i].property == Nnp[j])
      // test if number of property matchs number from input
      {
        if (Ncs[j])
        // test if any cross-section is assigned
        {
          fprintf(out, " %d %ld", Ncs[j], Dbcrs->crs[ic[i]].ridx[Ncsi[j]-1]);
          ncf = 1;
          break;
        }
      }
    }
    if (ncf == 0)
    // no cross section assigned
      fprintf(out, " 0");
    for (j = 0, nlsf = 0; j < nprop; j++)
    // loop over all properties from input
    {
      if (Top->N[i].property == Nnp[j])
      // test if number of property matchs number from input
      {
        if (Nlsv[j])
        // test if any local system is assigned
        {
          nlsf = 1;
          fprintf(out, " %ld", Nlsv[j][0].n);
          for (k = 0; k < Nlsv[j][0].n; k++)
          {
            for (l = 0; l < Nlsv[j][k].n; l++)
              fprintf(out, " % e", Nlsv[j][k](l));
          }
          fprintf(out, "\n");
          break;
        }
      }
    }
    if (nlsf == 0)
      fprintf(out, " 0\n");
  }
  fprintf(out, "\n");
  delete [] ic;
  return(0);
}



/**
 Function writes section with description of boundary conditions of nodes
  to the text file given by out.

 @param out - pointer to the opened text file, where the data will be written
 @param nprop - number of entries in the input property file

 @retval 0 - if succes
 @retval 1 - if nodes haven't assigned any boundary condition
*/
long wr_bocon(FILE *out, long nprop)
{
  long i, j, k, nbc, tcn;

  nbc = 0;
  tcn = 0;
  Nspd = 0;
  Ndpd = 0;
  Boccn = new long *[nprop];
  memset(Boccn, 0, nprop*sizeof(*Boccn));
  for (i = 0; i < Top->nodes; i++)
  {
    for (j = 0; j < nprop; j++)
    // loop over all properties from input
    {
      if (Top->N[i].property == Nnp[j])
      // test if number of property matchs number from input
      {
        if (Boc[j])
        // test if any boundary condition is prescribed
        {
          nbc++;
          if (Boccn[j] == NULL)
          {
            Boccn[j] = new long[Tndofn[i]];
            memset(Boccn[j], 0, Tndofn[i] * sizeof(long));
            for (k = 0; k < Tndofn[i]; k++)
            // loop for assigning bc code numbers
            // if direction has bc, increase cn
            {
              if ((Boc[j]->dir[k] >= 0) && (Boc[j]->con[k] != 0.0))
              {
                tcn++;
                Boccn[j][k] = tcn;
                if (Boc[j]->expr[k] == NULL)
                  Nspd++;
                else
                  Ndpd++;
              }
            }
          }
        }
        if (Nccn[j])
        // test if any common code numbers are prescribed
          nbc++;

      }
    }
  }
/*  if (nbc == 0)
  {
    fprintf(stderr, "\nError - problem hasn't assigned any boundary condition\n");
    return(1);
  }*/
  fprintf(out, "%ld\n", nbc);
  for (i = 0; i < Top->nodes; i++)
  {
    for (j = 0; j < nprop; j++)
    // loop over all properties from input
    {
      if (Top->N[i].property == Nnp[j])
      // test if number of property matchs number from input
      {
        if ((Boc[j]) || (Nccn[j]))
        // test if any boundary condition is prescribed or any common code numbers are prescribed
        {
          fprintf(out, "%ld", i+1);
          for (k = 0; k < Tndofn[i]; k++)
          // loop for writing bc
          // if direction has bc, write 0 else write 1
          {
            if (Boc[j] && Nccn[j])
            {
              if ((Boc[j]->dir[k] > 0) && (Nccn[j][k]))
              {
                fprintf(stderr, "\n\nError - node %ld has prescribed boundary condition and common code number\n", i+1);
                fprintf(stderr, " in the same direction\n.");
                fprintf(stderr, " function wr_bocon (input.cpp)\n.");
                return(2);
              }
              if ((Boc[j]->dir[k] < 0) && (Nccn[j][k] == 0))
              // no boundary condition nor common code number
              {
                fprintf(out, " 1");
                continue;
              }
            }
            if (Boc[j])
            {
              if (Boc[j]->dir[k] < 0)
              // no boundary condition
              {
                fprintf(out, " 1");
                continue;
              }
              if (Boc[j]->con[k] == 0.0)
              // rigid support
              {
                fprintf(out, " 0");
                continue;
              }
              if (Boc[j]->con[k] != 0.0)
              // prescribed displacement
              {
                fprintf(out, " %ld", -Boccn[j][k]);
                continue;
              }
            }
            if (Nccn[j])
            // common code number
            {
              if (Nccn[j][k])
              {
                fprintf(out, " %ld", Nccn[j][k]);
                continue;
              }
              if (Nccn[j][k] == 0)
              {
                fprintf(out, " 1");
                continue;
              }
            }
          }
          fprintf(out, "\n");
          break;
        }
      }
    }
  }
  fprintf(out, "\n");
  return(0);
}



/**
 Function writes section with description of elements to the text file given by out.
  @param out - pointer to the opened text file, where the data will be written
  @param nprop - number of entries in the input property file

  @retval 0 : on succes
  @retval 1 : if element hasn't assigned type
*/
long wr_elements(FILE *out, long nprop, descrip *d)
{
  long i, j, k, tip;
  crsectype *ce;
  long *ci, *ic;
  strastrestate ssst;
  long rindex;
  elemtype et;

  mattype **me;
  long **mi, **im, *nm;

  me   = new mattype* [Top->elements];
  memset(me, 0, sizeof(*me)*Top->elements);
  mi   = new long* [Top->elements];
  memset(mi, 0, sizeof(*mi)*Top->elements);
  im   = new long* [Top->elements];
  memset(im, 0, sizeof(*im)*Top->elements);
  nm   = new long  [Top->elements];

  Tet = new elemtype[Top->elements];
  memset(Tet, 0, sizeof(*Tet)*Top->elements);
  ce   = new crsectype [Top->elements];
  memset(ce, 0, sizeof(*ce)*Top->elements);
  ci   = new long  [Top->elements];
  memset(ci, 0, sizeof(*ci)*Top->elements);
  ic   = new long  [Top->elements];

  fprintf(out, "%ld\n", Top->elements);
  for (i = 0; i < Top->elements; i++)
  {
    for (j = 0, et = elemtype(-1); j < nprop; j++)
    // loop over all properties from input
    {
      // test if number of property matchs number from input
      if ((Top->E[i].property == Nep[j]) && (Sel[j] == Top->E[i].type))
      {
        if (Tel[j] >= 0)
        // assigning element type
        {
          et = Tel[j];
          ssst = Ssse[j];
          break;
        }
      }
    }
    if (et < 0)
    {
      fprintf(stderr, "\nError - element number %ld hasn't assigned type\n", i+1);
      for (k = 0; k < Top->elements; k++)
      {
        delete [] me[k]; delete [] mi[k]; delete [] im[k];
      }
      delete [] me; delete [] mi; delete [] nm;
      delete [] ce; delete [] ci; delete [] ic;
      return(1);
    }
    for (j = 0; j < nprop; j++)
    // loop over all properties from input
    {
      if (Ecs[j])
      // assigning cross-section to the element with number i
      {
        ce[i] = Ecs[j];
        ci[i] = Ecsi[j];
      }
    }
    if (ce[i] > 0)
    {
      if ((ic[i] = Dbcrs->search_crs(ce[i], ci[i]-1)) < 0)
      {
        fprintf(stderr, "\nError - cross-section with number %d and index %ld not found\n", ce[i], ci[i]);
        delete [] ce; delete [] ci; delete [] ic;
        return(1);
      }
      else
      {
        if (Dbcrs->crs[ic[i]].propu[ci[i]-1] == 0)
          Dbcrs->crs[ic[i]].npropu++; // indicator whether given cs type is used
        if (Ncu[ic[i]] == 0)
          Ncu[ic[i]]++; // indicator whether cs type from cs database is used in the given task
        Dbcrs->crs[ic[i]].propu[ci[i]-1]++;  // indicator whether given cs index is used
      }
    }
    for (j = 0; j < Dbcrs->numt; j++)
    {
      if (Dbcrs->crs[j].npropu)
      {
        rindex = 1;
        for (k = 0; k < Dbcrs->crs[j].nprop; k++)
        {
          if (Dbcrs->crs[j].propu[k])