!  Calculul aproximativ al fermelor cablu
!  Fisierul de date:
!  l,f1,f2,n    : semi-deschidere, sageti; nr. de panouri
!  A1,A2,A3   : arii - cablu sup., cablu inf., montanti
!  H          : forta de pretensionare - c. inf.
!  E,sigma    : modul de elast.; tensiune limita

program APROX_LIN

implicit real(k)      
real l,l2
character(260) fname, dummy*80

logical alta_incarcare
character(1) char1



K(E,A,l,f,H) = (E*A / (1+ .75*(l/f)**2) + H)*2/l**2
                  
print '(a\)',' Numele fisierului de date: '
read '(a)', fname
open(1, file=fname, status='old')

90 read(1,'(a)') dummy		! Antet
if(dummy(1:1) =='/' .and. dummy(2:2) =='/') then
goto 90
else
backspace 1
endif

read(1,*) l, f10, f20, n
read(1,*) a1, a2, a3
read(1,*) H20
read(1,*) E,sigma

l2 = l*l
a = 2*l/n
h10 = h20*f20/f10
p0 = h20*f20*2/l2

K1 = K(E,A1,l,f10,H10)
K2 = K(E,A2,l,f20,H20)
write(1,'(/, a)') 'Rigiditati:'
write(1,*) ' K1 = ',K1
write(1,*) ' K2 = ',K2      
write(1,*)

! incarcarea initiala in montanti = suctiune de slabire cablu sup.
print*,'  Incarcare echivalenta din pretensionare:', p0
write(1,*) '  Incarcare echivalenta din pretensionare:', p0

! incarcarea de slabire montanti: dp=p0
q_slabire=p0*(K1+K2)/K2
print*,' *Incarcare de slabire montanti          :', q_slabire
write(1,*) ' *Incarcare de slabire montanti          :', q_slabire

! incarcarea de slabire cablu superior:
q_slab_cs =-p0*(K1+K2)/K1
print*,' *Incarcare de slabire cablu superior    :', q_slab_cs
write(1,*) ' *Incarcare de slabire cablu superior    :', q_slab_cs
      

iInc =0
100 iInc =iInc+1
icaz=0
print'(/, a, /, a)',' Introdu incarcarea: q_min, q_max, q_pas:', ' (q, 0, 0 = pt. o singura incarcare)'
read*, q_min, q_max, q_pas
if(q_max ==0. .and. q_pas ==0.) then
q_max =q_min
q_pas =1.
endif
!  Write results header:
write(1, '(/, a, I3, a)') 'INCARCAREA nr.', iInc, ':'

write(1,10)' Nr.crt.','Incarcare', 'Sageata "df"', 'H1',       'H2',      'T-montanti'
10 format(/, t1,a,     t12,a,       t24,a,          t39,a,      t52,a,     t65,a)
11 format(/, I5,       t11,f10.3,   t24,f10.6,      t36,f11.4,  t49,f11.4, t62,f10.4, \)

DO q=q_min, q_max, q_pas
icaz = icaz+1

dp = q*K2/(K1+K2)
df = q / (K1+K2)
      
p  = p0  -dp
f1 = f10 +df
f2 = f20 -df
      
H2 = p*l2/2/f2
H1 = (p+q)*l2/2/f1

! Varianta:
!p$q =p0 +K1*df
!H1 = p$q*l2/2/f1

if(q ==0.)then                                
icaz_cur=0
icaz=icaz-1
else
icaz_cur=icaz
endif

!  write results:
write(1,11) icaz_cur, q, df, H1, H2, p*a
if(H1/A1 >=sigma) then
	print*,' * Cablu sup. rupt'
	write(1,*) '   * Cablu sup. rupt'
elseif(H2/A2.ge.sigma) then
	print*,' * Cablu inf. rupt'
    write(1,*) '   * Cablu sup. rupt'
elseif(p.lt.0) then
	print*,' * Montanti slabiti'
    write(1,*) '   * Montanti slabiti'
elseif(p*a/A3 >=sigma)then
	print*,' * Montanti rupti'
    write(1,*) '   * Montanti rupti'
endif
write(1,*)
ENDDO

print '(/, a, \)', 'Alta incarcare ? (Y, y/ char): '
read (*, '(a1)') char1

alta_incarcare =char1 =='Y' .or.char1 =='y'

if(alta_incarcare) then
goto 100
endif

500   close (1)

write(*,'(/, a, a, /)'),' * Fisier de rezultate :', trim(fname)

END

