c==============================================================================
c=== PROGRAMA DE CALCULO POR ELEMENTOS FINITOS PARA PROBLEMAS ELIPTICOS 1D ====
c==============================================================================

c______________________________________Programa MEF1D__________________________
c
c AUTORES: F. Navarrina, I. Colominas, H. Gomez, M. Casteleiro
c                                                                 _
c        Dep. de Metodos Matematicos y de Rep., ETSICCP de La Coruna, UDC.
c
c FECHA ULTIMA REVISION: 20070524
c
c Copyright: Este programa es propiedad de su autor.
c            No puede utilizarse con fines comerciales sin conocimiento
c            y autorizacion de su autor.
c            Puede utilizarse con fines educativos sin ninguna restriccion.
c
c------------------------------------------------------------------------------
c===============================================================================
c===============================================================================
c
c   MEF1D - Metodo de Elementos Finitos 1D:
c                                           (0) Elasticidad Lineal 1D
c
c===============================================================================
c===============================================================================

      implicit integer*4 (i-n)
d     implicit real*8 (a-h,o-z)

      parameter(mxdime=    1,mxdofn=    1,mxstre=    1,mxprop=    5,
     .          mxnode=    3,mxgaue=    3,mxmats=  100,mxvfix=  200,
     .          mxload=    5,mxcase=    5,
     .          mxelem=10000,mxpoin=20001,mxperf=  200)

      parameter(mxevab=mxdofn*mxnode,mxtotv=mxdofn*mxpoin,
     .          mxtotc=mxcase*mxtotv,mxrigi=mxperf*mxtotv)

      common /PAR/ntype,nelem,nnode,npoin,ndime,ndofn,nstre,nmats,nprop,
     .            ngaue,nload,ncase,nvfix,nevab,ntotv,ntotc,nrigi
   
      character*40 archivo
      character    titul*80,tcase(mxcase)*40
      character*1  siono

c_____Definicion del problema
      dimension lnods(mxnode*mxelem),matno(mxelem),props(mxprop*mxmats),
     .          nofix(2*mxvfix),presc(mxvfix),ngfix(mxvfix),
     .          coord(mxdime*mxpoin)

c_____Nivel de elemento
      dimension posge(mxdime*mxgaue),weige(mxgaue),
     .          shape(mxnode*mxgaue),derie(mxdime*mxnode*mxgaue),
     .          elcod(mxdime*mxnode),gcord(mxdime*mxgaue),
     .          xjace(mxdime*mxdime*mxgaue),djacb(mxgaue),
     .          xjaci(mxdime*mxdime*mxgaue),cartd(mxdime*mxnode*mxgaue),
     .          dmatx(mxstre*mxstre),bmatx(mxstre*mxevab*mxgaue),
     .          smatx(mxstre*mxevab*mxgaue),dvolu(mxgaue),
     .          estif(mxevab*mxevab),eload(mxevab),eldis(mxdofn*mxnode),
     .          strsg(mxstre*mxgaue*mxcase)

c_____Tratamiento de cargas
      dimension fload(mxload*mxcase),kload(mxload),nilod(mxload),
     .          rilod(3,mxload),pforc(mxdofn),bforc(mxdofn*mxgaue)

c_____Reacciones
      dimension treac(mxvfix*mxcase),totre(mxdofn)

c_____Matriz de rigidez y vector de fuerzas=vector de desplazamientos
      dimension lpntn(mxpoin),lpamc(mxtotv)
      dimension gload(mxtotc),rgmat(mxrigi)

c===============================================================================
c===============================================================================

c.....COMIENZO 
      call SRV$Get_CpuTime(itime_0)

C.....UNIDADES LOGICAS DE E/S
      logd=10
      logr=20
      loge= 6

C.....ARCHIVOS DE DATOS Y RESULTADOS
    1 write(6,'(a,$)') ' > Nombre Generico de Archivos : '
      read(5,'(q,a)') lon,archivo(1:lon)
      open(unit=logd,file=archivo(1:lon)//'.dat',status='old',err=1)
      open(unit=logr,file=archivo(1:lon)//'.res',status='new',err=1)
 
C.....DEFINICION DE MALLA Y PARAMETROS DE CONTROL
      read(logd,'()')
      read(logd,'(i6)') ntype
      read(logd,'()')
      read(logd,'(1x,a)')  titul
      read(logd,'()')
      read(logd,'(12i6)') npoin,nvfix,nelem,nnode,nmats,ngaue
      if (npoin.gt.mxpoin) call SRV$Error('MEF1D: '//
     .  'Memoria insuficiente (NPOIN > MXPOIN)')
      if (nvfix.gt.mxvfix) call SRV$Error('MEF1D: '//
     .  'Memoria insuficiente (NVFIX > MXVFIX)')
      if (nelem.gt.mxelem) call SRV$Error('MEF1D: '//
     .  'Memoria insuficiente (NELEM > MXELEM)')
      if (nnode.gt.mxnode) call SRV$Error('MEF1D: '//
     .  'Memoria insuficiente (NNODE > MXNODE)')
      if (nmats.gt.mxmats) call SRV$Error('MEF1D: '//
     .  'Memoria insuficiente (NMATS > MXMATS)')
      if (ngaue.gt.mxgaue) call SRV$Error('MEF1D: '//
     .  'Memoria insuficiente (NGAUE > MXGAUE)')

      call CAL$Ini_Par(ntype,ndime,npoin,ndofn,nnode,nstre,ngaue,
     .                 nprop,ntloa,nevab,ntotv,posge,weige,shape,derie)

      call PAR$Def_Malla(logd,ndime,npoin,coord,lpntn,nnode,nelem,lnods,
     .                   matno,nprop,nmats,props,nvfix,nofix,presc)

C.....DEFINICION DE TIPOS E HIPOTESIS DE CARGA
      read(logd,'()')
      read(logd,'(2i6)') ncase,nload
      if (ncase.gt.mxcase) call SRV$Error('MEF1D: '//
     .  'Memoria insuficiente (NCASE > MXCASE)')
      if (nload.gt.mxload) call SRV$Error('MEF1D: '//
     .  'Memoria insuficiente (NLOAD > MXLOAD)')

      call PAR$Def_Carga(logd,ncase,nload,tcase,fload,kload,nilod,rilod)

C.....INICIALIZACION DE MEMORIA
      call CAL$Ini_AMC(nelem,nnode,npoin,ndofn,ntotv,lnods,lpamc,
     .                 lpntn,nrigi)
      if (nrigi.gt.mxrigi) call SRV$Error('MEF1D: '//
     .  'Memoria insuficiente (NRIGI > MXRIGI)')

      call SRV$Ini_Vector(nrigi,rgmat)  !matriz de rigidez
      ntotc=ncase*ntotv
      call SRV$Ini_Vector(ntotc,gload)  !vector de fuerzas nodales

C.....ESCRITURA DE DATOS
      call PAR$Rotular(logr,'    BLOQUE 1  --->  DATOS        ')

      call PAR$Esc_Datos(logr,ntype,titul,ndime,ndofn,nstre,ntloa,
     .     npoin,nnode,nelem,nprop,nmats,nload,ncase,ngaue,nvfix,
     .     nevab,ntotv,ntotc,nrigi,coord,lpntn,lnods,matno,props,
     .     nofix,presc,kload,nilod,rilod,fload,tcase)
 
C.....LECTURA, ENSAMBLAJE Y ESCRITURA DE CARGAS PUNTUALES
      call PAR$Def_FPun(logd,logr,ndofn,npoin,ntotv,nload,ncase,
     .     kload,nilod,pforc,fload,lpntn,gload)

      close(unit=logd)

C.....GENERACION DE MATRIZ DE RIGIDEZ Y VECTORES DE FUERZAS
      call SRV$Get_CpuTime(itime_1)

c.....Bucle sobre elementos
      do ielem=1,nelem

c.......Generacion, integracion y ensamblaje de matriz de rigidez de elementos
        call CAL$Proc_Rig(2,ntype,ielem,nelem,nnode,npoin,ndime,
     .       ndofn,nstre,nmats,nprop,ngaue,nevab,ntotv,nrigi,lnods,
     .       coord,elcod,matno,props,weige,shape,derie,gcord,dvolu,
     .       xjace,djacb,xjaci,cartd,dmatx,bmatx,smatx,estif,lpamc,
     .       lpntn,rgmat)

c.......Generacion, integracion y ensamblaje de vectores de fuerzas volumetricas
        do iload=1,nload
          call CAL$Proc_FVol(ntype,ielem,iload,nelem,nnode,npoin,
     .         ndime,ndofn,nmats,nprop,ngaue,ntotv,nload,ncase,lnods,
     .         matno,props,shape,dvolu,bforc,gcord,kload,rilod,fload,
     .         eload,lpntn,gload)
        enddo

      enddo

c.....imposicion de condiciones de vinculacion
      call CAL$Imp_Vincul(nvfix,ndofn,npoin,nofix,ngfix,lpntn)

C.....SOLUCION DEL SISTEMA DE ECUACIONES (Metodo de Cholesky)
      call SRV$Get_CpuTime(itime_2)
      call SLE$Solver_LDLt_CP(0,
     .                        ntotv,lpamc,rgmat,
     .                        nvfix,ngfix,presc,
     .                        ncase,gload,treac)

C.....ESCRITURA DE RESULTADOS
      call SRV$Get_CpuTime(itime_3)

c.....Desplazamientos y reacciones
      do icase=1,ncase
        call PAR$Rotular(logr,'   BLOQUE 2  --->  DESPL./REACC. ')


        write(logr,2100) icase,tcase(icase)
 2100   format(//
     .    /1x,130('=')
     .    /1x,'HIPOTESIS DE CARGA NUMERO ',i5,'    --------->     ',a40
     .    /1x,130('='))

        call PAR$Esc_DR(logr,icase,ndofn,npoin,nvfix,nofix,lpntn,
     .                  ncase,gload,treac,totre)
      enddo

c.....Tensiones en puntos de integracion
      call PAR$Rotular(logr,'    BLOQUE 3  --->  TENSIONES    ')

c.....Bucle sobre elementos
      do ielem=1,nelem

        write(logr,3100) ielem
 3100   format(//
     .    /1x,130('=')
     .    /1x,'ELEMENTO NUMERO           ',i8
     .    /1x,130('='))

        call CAL$Proc_Rig(0,ntype,ielem,nelem,nnode,npoin,ndime,
     .       ndofn,nstre,nmats,nprop,ngaue,nevab,ntotv,nrigi,lnods,
     .       coord,elcod,matno,props,weige,shape,derie,gcord,dvolu,
     .       xjace,djacb,xjaci,cartd,dmatx,bmatx,smatx,estif,lpamc,
     .       lpntn,rgmat)

        call CAL$Proc_Tens(ielem,nelem,ncase,nnode,ndofn,npoin,
     .       nstre,nevab,ngaue,lnods,lpntn,gload,eldis,smatx,strsg)

        call PAR$Esc_Tens(logr,nstre,ngaue,ncase,strsg)

      enddo

c.....RESUMEN DE TIEMPO DE COMPUTACION UTILIZADO
      call SRV$Get_CpuTime(itime_4)

      itime_tot=itime_4-itime_0
      ttime_tot=dfloat(itime_tot)
      tsegs=ttime_tot/.1d+03

      call PAR$Rotular(logr,'  BLOQUE 4  --->  TIEMPO DE CPU  ')

      write(logr,4100) tsegs
 4100 format(//
     .   /1x,15x,'TIEMPO TOTAL DE CPU ..........',f15.2,' Segundos')

      if (itime_tot.ne.0) then

        ptime_ini=dfloat(itime_1-itime_0)/ttime_tot*.1d+03
        ptime_gen=dfloat(itime_2-itime_1)/ttime_tot*.1d+03
        ptime_sol=dfloat(itime_3-itime_2)/ttime_tot*.1d+03
        ptime_esc=dfloat(itime_4-itime_3)/ttime_tot*.1d+03

        write(logr,4110) ptime_ini,ptime_gen,ptime_sol,ptime_esc
 4110   format(/1x,15x,'Inicializacion ...............',f15.2,' %'
     .         /1x,15x,'Generacion e Integracion .....',f15.2,' %'
     .         /1x,15x,'Solucion de Sistemas .........',f15.2,' %'
     .         /1x,15x,'Escritura de Resultados ......',f15.2,' %')

      endif

      close (unit=logr)

      end

!===============================================================================
!===============================================================================
!                               INPUT/OUTPUT
!===============================================================================
!===============================================================================

c__________________________________________________________________PAR$Def_Malla
c     Lectura de Propiedades, Coordenadas, Conectividad y Coacciones

      subroutine PAR$Def_Malla(logd,ndime,npoin,coord,lpntn,nnode,nelem,
     .           lnods,matno,nprop,nmats,props,nvfix,nofix,presc)

      implicit integer*4 (i-n)
d     implicit real*8 (a-h,o-z)
   
      dimension coord(ndime,npoin),lpntn(npoin),lnods(nnode,nelem),
     .          matno(nelem),props(nprop,nmats),nofix(2,nvfix),
     .          presc(nvfix)

c.....Propiedades de material
      read(logd,'()')
      do imats=1,nmats
        read(logd,101) jmats,(props(iprop,jmats),iprop=1,nprop)
  101   format(i6,5e20.9)
      enddo
 
c.....Coordenadas nodales
      read(logd,'()')
      do ipoin=1,npoin
        read(logd,111) jpoin,lpntn(jpoin),
     .                    (coord(idime,jpoin),idime=1,ndime)
  111 format(2i6,3e20.9)
      enddo

c.....Conectividad y caracteristicas del elemento
      read(logd,'()')
      do ielem=1,nelem
        read(logd,121) jelem,knode,matno(jelem),
     .                    (lnods(inode,jelem),inode=1,nnode)
  121   format(3i6,10i6)
      enddo

c.....Coacciones
      read(logd,'()')
      do ivfix=1,nvfix
        read(logd,131) jvfix,nofix(1,jvfix),nofix(2,jvfix),presc(jvfix)
  131   format(3i6,e20.9)
      enddo

      return
      end

c__________________________________________________________________PAR$Def_Carga
c     Lectura de Estados de carga

      subroutine PAR$Def_Carga(logd,ncase,nload,tcase,fload,kload,nilod,
     .           rilod)

      implicit integer*4 (i-n)
d     implicit real*8 (a-h,o-z)

      character*(*) tcase(ncase)
      dimension fload(ncase,nload),kload(nload),nilod(nload),
     .          rilod(3,nload)

c.....Hipotesis de carga
      read(logd,'()')
      do icase=1,ncase
        read(logd,'(i6,1x,a)') jcase,tcase(jcase)
        read(logd,'(5e20.9)') (fload(jcase,iload),iload=1,nload)
      enddo

c.....Tipos de carga
      read(logd,'()')
      do iload=1,nload
        read(logd,'(3i6,3e20.9)') jload,kload(jload),
     .                            nilod(jload),(rilod(i,jload),i=1,3)
      enddo
        
      return
      end

c____________________________________________________________________PAR$Rotular
c     Escribe Rotulo de Bloque de Resultados

      subroutine PAR$Rotular(logr,rotul)
      character*(*) rotul

      write(logr,1000) rotul
 1000 format(
     .   '1',131('*')
     .   /1x,131('*')
     .   /1x, 45('*'),4x,'                                 ',4x,45('*')
     .   /1x, 45('*'),4x,'MEF/1D - PROBLEMA ELASTICO LINEAL',4x,45('*')
     .   /1x, 45('*'),4x,'                                 ',4x,45('*')
     .   /1x, 45('*'),4x,              a33                  ,4x,45('*')
     .   /1x, 45('*'),4x,'                                 ',4x,45('*')
     .   /1x,131('*')
     .   /1x,131('*'))

      return
      end

c__________________________________________________________________PAR$Esc_Datos
c     Escritura de datos

      subroutine PAR$Esc_Datos(logr,ntype,titul,ndime,ndofn,nstre,ntloa,
     .           npoin,nnode,nelem,nprop,nmats,nload,ncase,ngaue,nvfix,
     .           nevab,ntotv,ntotc,nrigi,coord,lpntn,lnods,matno,props,
     .           nofix,presc,kload,nilod,rilod,fload,tcase)

      implicit integer*4 (i-n)
d     implicit real*8 (a-h,o-z)

      character*(*) titul

      dimension coord(ndime,npoin),lpntn(npoin),lnods(nnode,nelem),
     .          matno(nelem),props(nprop,nmats),nofix(2,nvfix),
     .          presc(nvfix),kload(nload),nilod(nload),rilod(3,nload),
     .          fload(ncase,nload)

      character*(*) tcase(ncase)

c.....datos de control
      write(logr,1110) titul
 1110 format(//
     .  /1x,130('=')
     .  /1x,130('=')
     .  /1x,'TITULO DEL PROBLEMA                --------->     ',a80
     .  /1x,130('=')
     .  /1x,130('='))

      write(logr,1200)
 1200 format(//
     .  /1x,'DATOS DE CONTROL              ',100('.')
     .  /1x,130('-'))
      write(logr,1210) ntype,nelem,nnode,npoin,ndime,ndofn,
     .                 nstre,nmats,nprop,ngaue,nload,ncase,
     .                 nvfix,nevab,ntotv,ntotc,nrigi
 1210 format(/1x,15x,'Ntype =',i8,5x,'Nelem =',i8,5x,'Nnode =',i8,
     .            5x,'Npoin =',i8,5x,'Ndime =',i8,5x,'Ndofn =',i8,
     .       /1x,15x,'Nstre =',i8,5x,'Nmats =',i8,5x,'Nprop =',i8,
     .            5x,'Ngaue =',i8,5x,'Nload =',i8,5x,'Ncase =',i8,
     .       /1x,15x,'Nvfix =',i8,5x,'Nevab =',i8,5x,'Ntotv =',i8,
     .            5x,'Ntotc =',i8,5x,'Nrigi =',i8)

      write(logr,1300)
 1300 format(//
     .  /1x,'MATERIALES Y PROPIEDADES      ',100('.')
     .  /1x,130('-'))
      write(logr,1310)
 1310 format(/1x,15x,'Imats',15x,'Young',15x,'Poiss',15x,'Thick',
     .           15x,'Dense',15x,'Alpha'/)
      do imats=1,nmats
        write(logr,1320) imats,(props(iprop,imats),iprop=1,nprop)
 1320   format(1x,i20,5e20.9)
      enddo

      write(logr,1400)
 1400 format(//
     .  /1x,'COORDENADAS NODALES           ',100('.')
     .  /1x,130('-'))
      write(logr,1410)
 1410 format(/1x,15x,'Ipoin',5x,'Lpntn',15x,
     .               '(Coord(idime,ipoin),idime=1,ndime)'/)
      do ipoin=1,npoin
        write(logr,1420) ipoin,lpntn(ipoin),
     .                 (coord(idime,ipoin),idime=1,ndime)
 1420   format(1x,10x,2i10,10x,3e20.9)
      enddo

      write(logr,1500)
 1500 format(//
     .  /1x,'CONECTIVIDAD, TIPO DE ELEMENTO',100('.')
     .  /1x,130('-'))
      write(logr,1510)
 1510 format(/1x,15x,'Ielem',5x,'Nnode',5x,'Matno',15x,
     .               '(Lnods(inode,ielem),inode=1,nnode)'/)

      do ielem=1,nelem
        write(logr,1520) ielem,nnode,matno(ielem),
     .                  (lnods(inode,ielem),inode=1,nnode)
 1520     format(1x,10x,3i10,10x,9i8)
      enddo

      write(logr,1600)
 1600 format(//
     .  /1x,'COACCIONES                    ',100('.')
     .  /1x,130('-'))
      write(logr,1610)
 1610 format(/1x,15x,'Ivfix',5x,'Ipoin',5x,'Idofn',15x,'Presc'/)
      do ivfix=1,nvfix
        write(logr,1620) ivfix,nofix(1,ivfix),nofix(2,ivfix),
     .                   presc(ivfix)
 1620   format(1x,10x,3i10,10x,e20.9)
      enddo

      write(logr,1700)
 1700 format(//
     .  /1x,'HIPOTESIS DE CARGA            ',100('.')
     .  /1x,130('-'))
      do icase=1,ncase
        write(logr,1702) icase,tcase(icase)
 1702   format(/1x,10x,'Hipotesis Numero ',i3,
     .             5x,'---->',5x,a40,45('.')
     .         /1x,10x,20('-'),15(' '),85('-'))
        write(logr,1710)
 1710   format(/1x,15x,'Iload (carga elemental)       ',
     .                 'Fload (coeficiente de seguridad)'/)
        do iload=1,nload
          write(logr,1720) iload,fload(icase,iload)
 1720     format(1x,10x,i10,20x,e20.9)
        enddo
      enddo

      write(logr,1800)
 1800 format(//
     .  /1x,'CARGAS ELEMENTALES            ',100('.')
     .  /1x,130('-'))
      write(logr,1810)
 1810 format(/1x,15x,'Iload',5x,'Kload',10x,'Descripcion'/)
      do iload=1,nload
        if (kload(iload).eq.1) then
          write(logr,1820) iload,kload(iload),
     .                     'PESO PROPIO                   ',
     .                     'aceleracion gravitatoria    = ',
     .                     rilod(1,iload)
        else if (kload(iload).eq.2) then
          w=rilod(2,iload)*rilod(2,iload)+rilod(3,iload)*rilod(3,iload)
          w=dsqrt(w)
          write(logr,1820) iload,kload(iload),
     .                     'FUERZA CENTRIFUGA             ',
     .                     'velocidad de rotacion       = ',
     .                      w
        else if (kload(iload).eq.3) then
          write(logr,1822) iload,kload(iload),
     .                     'FUERZAS PUNTUALES             ',
     .                     'numero de fuerzas aplicadas = ',
     .                      nilod(iload)
        endif
 1820   format(1x,10x,2i10,10x,a,a,e20.9)
 1822   format(1x,10x,2i10,10x,a,a,i10)
      enddo

      return
      end

c___________________________________________________________________PAR$Def_FPun
c     Lectura, ensamblaje y escritura de cargas puntuales

      subroutine PAR$Def_FPun(logd,logr,ndofn,npoin,ntotv,nload,ncase,
     .           kload,nilod,pforc,fload,lpntn,gload)

      implicit integer*4 (i-n)
d     implicit real*8 (a-h,o-z)

      dimension kload(nload),nilod(nload),pforc(ndofn)

c.....bucle sobre cargas tipo
      do iload=1,nload

        if (kload(iload).eq.3) then                              ! Carga Puntual

          read (logd,'()')
          read (logd,'(i6)') jload

          write(logr,1902) jload
 1902     format(//
     .      /1x,10x,'Sistema de Cargas Puntuales Numero',i6
     .      /1x,10x,40('-'))
          write(logr,1910)
 1910     format(/1x,15x,'Iilod',5x,'Ipoin',15x,
     .               '(Pforc(idofn),idofn=1,ndofn)'/)

          do iilod=1,nilod(jload)
            read (logd,1921) jilod,jpoin,(pforc(idofn),idofn=1,ndofn)
 1921       format(6x,2i6,3e20.9)
            write(logr,1920) jilod,jpoin,(pforc(idofn),idofn=1,ndofn)
 1920       format(1x,10x,2i10,10x,3e20.9)

            call CAL$Proc_FPun(jpoin,jload,ndofn,npoin,ntotv,
     .           nload,ncase,pforc,fload,lpntn,gload)
          enddo

        endif

      enddo    ! fin de bucle sobre estados de carga

      return
      end

c_____________________________________________________________________PAR$Esc_DR
c     Escritura de resultados (desplazamientos y reacciones)

      subroutine PAR$Esc_DR(logr,icase,ndofn,npoin,nvfix,nofix,lpntn,
     .                      ncase,xdisp,treac,totre)

      implicit integer*4 (i-n)
d     implicit real*8 (a-h,o-z)

      dimension lpntn(npoin),nofix(2,nvfix),xdisp(ndofn,npoin,ncase),
     .          treac(nvfix,ncase),totre(ndofn)

c.....Desplazamientos
      write(logr,2200)
 2200 format(//
     .  /1x,'DESPLAZAMIENTOS NODALES       ',100('.')
     .  /1x,130('-'))
      write(logr,2210)
 2210 format(/1x,15x,'Ipoin',5x,'Lpntn',15x,
     .       '(Xdisp(idofn,Lpntn(ipoin),icase),idofn=1,ndofn)'/)
      do ipoin=1,npoin
        ipoir=lpntn(ipoin)
        write(logr,2220) ipoin,ipoir,
     .       (xdisp(idofn,ipoir,icase),idofn=1,ndofn)
 2220   format(1x,10x,2i10,10x,3e20.9)
      enddo

c.....Inicializacion de resultante de reacciones
      call srv$ini_vector(ndofn,totre)

c.....Reacciones
      write(logr,2300)
 2300 format(//
     .  /1x,'REACCIONES                    ',100('.')
     .  /1x,130('-'))
      write(logr,2310)
 2310 format(/1x,15x,'Ivfix',5x,'Ipoin',5x,'Idofn',15x,
     .       'Treac(ivfix,icase)'/)
      do ivfix=1,nvfix
        write(logr,2320) ivfix,nofix(1,ivfix),nofix(2,ivfix),
     .                   treac(ivfix,icase)
 2320   format(1x,10x,3i10,10x,e20.9)
        totre(nofix(2,ivfix))=totre(nofix(2,ivfix))+treac(ivfix,icase)
      enddo

c.....Resultante de Reacciones
      write(logr,2400)
 2400 format(//
     .  /1x,'RESULTANTE DE LAS REACCIONES  ',100('.')
     .  /1x,130('-'))
      write(logr,2410)
 2410 format(/1x,15x,'Idofn',15x,'Reaccion(idofn)'/)
      do idofn=1,ndofn
        write(logr,2420) idofn,totre(idofn)
 2420   format(1x,10x,i10,10x,e20.9)
      enddo

      return
      end

c___________________________________________________________________PAR$Esc_Tens
c     Escritura de Tensiones en puntos de integracion

      subroutine PAR$Esc_Tens(logr,nstre,ngaue,ncase,strsg)

      implicit integer*4 (i-n)
d     implicit real*8 (a-h,o-z)

      dimension strsg(nstre,ngaue,ncase)

      do icase=1,ncase

        write(logr,3200) icase
 3200   format(//
     .    /1x,10x,'Tensiones en la Hipotesis de Carga Numero : ',i6
     .    /1x,10x,50('-'))
        write(logr,3210)
 3210   format(/1x,15x,'Igaue',15x,
     .                 '(Strsg(istre,igaue,icase),istre=1,nstre)'/)
        do igaue=1,ngaue
          write(logr,3220) igaue,
     .                     (strsg(istre,igaue,icase),istre=1,nstre)
 3220     format(1x,10x,i10,10x,3e20.9)
        enddo

      enddo

      return
      end
