MODULE ModMES
contains
SUBROUTINE rozw(A,B)
IMPLICIT NONE
! argumenty
REAL, INTENT(IN OUT) :: A(:,:) ! macierz
REAL, INTENT(IN OUT) :: B(:) ! wek.prawych stron i rozw.
INTEGER:: i,j,n ! n -rzad macierzy A
LOGICAL:: ok
REAL :: element, Wiersz(SIZE(B,1)),pivtol
n=SIZE(B,1)
ok=SIZE(A,1)== n .AND. SIZE(A,2)== n
IF (.NOT.ok) THEN
WRITE(*,*)' ****** Niewlasciwe wymiary macierzy ****** '
STOP
END IF
pivtol=MAXVAL(A)*1e-10
DO j=1,n
DO i=1,j-1
A(i+1:n,j)=A(i+1:n,j)-A(i,j)*A(i+1:n,i)
END DO
i=MAXVAL(MAXLOC(ABS(A(j:n,j))))+j-1
IF (ABS(A(i,j)) < pivtol) THEN
ok=.FALSE.
RETURN
END IF
IF (i.NE.j) THEN
Wiersz=A(j,:)
A(j,:)=A(i,:)
A(i,:)=Wiersz
element=B(j)
B(j)=B(i)
B(i)=element
END IF
A(j+1:n,j)=A(j+1:n,j)/A(j,j)
END DO
DO i=1,n-1
B(i+1:n)=B(i+1:n)-B(i)*A(i+1:n,i)
END DO
DO j=n,1,-1
B(j)=B(j)/A(j,j)
B(1:j-1)=B(1:j-1)-B(j)*A(1:j-1,j)
END DO
RETURN
END SUBROUTINE rozw
!-----------------------------------------------------------
SUBROUTINE agre(w1,w2,Ke,Kg)
! agreg. mac. Ke do Kg blokami odpowiadajacymi wezlom w1 i w2
IMPLICIT NONE
REAL, INTENT(IN) :: Ke(:,:)
REAL, INTENT(IN OUT):: Kg(:,:)
INTEGER, INTENT(IN) :: w1,w2
INTEGER:: i,j,lsw,nwg,nkg,nwe,nke,W(2)
lsw=UBOUND(Ke,1)/2
W(1)=w1
W(2)=w2
DO i=1,2
DO j=1,2
nwg= (W(i)-1)*lsw+1
nkg= (W(j)-1)*lsw+1
nwe= (i-1)*lsw+1
nke= (j-1)*lsw+1
Kg(nwg:nwg+lsw-1,nkg:nkg+lsw-1)= Kg(nwg:nwg+lsw-1,nkg:nkg+lsw-1)+&
Ke(nwe:nwe+lsw-1,nke:nke+lsw-1)
END DO
END DO
RETURN
END SUBROUTINE agre
!---------------------------------------
END MODULE ModMES
PROGRAM belka
USE ModMes
IMPLICIT NONE
REAL::Ke(4,4),ej,kg(8,8),l(3),Fg(8),ue(4),Se(4)
INTEGER::e
KG=0.
l=(/2.0,3.0,2.5/)
FG=(/0.0,0.0,52.5,0.0,0.0,0.0,70.0,0.0/)
ej=6000.0
Do e=1,3
Ke(1,:)=(/(12.0*ej)/l(e)**3,(6.0*ej)/l(e)**2,(-12.0*ej)/l(e)**3,(6.0*ej)/l(e)**2/)
Ke(2,:)=(/(6.0*ej)/l(e)**2,(4.0*ej)/l(e),(-6.0*ej)/l(e)**2,(2.0*ej)/l(e)/)
Ke(3,:)=(/(-12.0*ej)/l(e)**3,(-6.0*ej)/l(e)**2,(12.0*ej)/l(e)**3,(-6.0*ej)/l(e)**2/)
Ke(4,:)=(/(6.0*ej)/l(e)**2,(2.0*ej)/l(e),(-6.0*ej)/l(e)**2,(4.0*ej)/l(e)/)
CALL agre(e,e+1,ke,kg)
END DO
kg(1,:)=0.
kg(:,1)=0.
kg(1,1)=1.
kg(2,:)=0.
kg(:,2)=0.
kg(2,2)=1.
kg(5,:)=0.
kg(:,5)=0.
kg(5,5)=1.
fg(1)=0.
fg(2)=0.
fg(5)=0.
CALL rozw(kg,fg)
OPEN (UNIT=1,FILE='WYNIKI.TXT',STATUS='replace',ACTION='WRITE')
WRITE(1,*)'PRZEMIESZCZENIA I KAT OBROTU/WARTOŚCI SIL WEWNETRZNYCH'
WRITE(1,*)FG
DO E=1,3
Ke(1,:)=(/(12.0*ej)/l(e)**3,(6.0*ej)/l(e)**2,(-12.0*ej)/l(e)**3,(6.0*ej)/l(e)**2/)
Ke(2,:)=(/(6.0*ej)/l(e)**2,(4.0*ej)/l(e),(-6.0*ej)/l(e)**2,(2.0*ej)/l(e)/)
Ke(3,:)=(/(-12.0*ej)/l(e)**3,(-6.0*ej)/l(e)**2,(12.0*ej)/l(e)**3,(-6.0*ej)/l(e)**2/)
Ke(4,:)=(/(6.0*ej)/l(e)**2,(2.0*ej)/l(e),(-6.0*ej)/l(e)**2,(4.0*ej)/l(e)/)
UE=FG((E-1)*2+1:(E-1)*2+4)
SE=MATMUL(KE,UE)
WRITE(1,*)E,SE
END DO
STOP
END PROGRAM BELKA