##Naide2_5.m
#  vt lk 55
close all 
clear
## 
#dbwhere(251)
#
##=========================================================================
## This Program is written by Andres Lahe,   2015-03-30
##                    e-mail: andres.lahe@ttu.ee
## LAST MODIFIED: Andres Lahe,    2015-05-05
## Copyright (c)  2015  by Tallinn University of Technology 
##                Department of Mechanics
##                http://www.ttu.ee/%
##
## This program is free software; you can redistribute it and/or
## modify it under the terms of the GNU General Public License
## as published by the Free Software Foundation; either version 2
## of the License, or (at your option) any later version.

## This program is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
## GNU General Public License for more details.

## You should have received a copy of the GNU General Public License
## along with this program; if not, write to the Free Software
## Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
## http://www.gnu.org/licenses/old-licenses/gpl-2.0.html
##=========================================================================
#
# eps2=0.0000000001
 eps1=1.0*10e-11;
 eps2=0.00001 # 4 nulli pärast koma, eristab sortimisel  199.99999
 global spA
 global X
 global AP
NNK=8 
#--------- Define a matrix spA to be sparse -----------
disp('-------- Sparse matrix instantiation --------  ')
disp('  spA=sparse(NNK,NNK)  ')
spA=sparse(NNK,NNK) # nullistan võrrandisüsteemi kordajate maatriksi
#--------- Define a vector B ------------
disp('-------- Right-hand side of the equations (RHS). --------  ')
disp('  Bvb=zeros(NNK,1);  ')
Bvb=zeros(NNK,1);   # nullistan võrrandisüsteemi vabaliikmed
disp('   ')
#######
#baasi0=1.0   # baasjäikus
l=400  # cm   # varda pikkus     
#x=0.0 # cm   # koordinaat x
a=200 ; # cm   # koondmomendi algus
Mx=10000.0 # Ncm    # koondmoment 
#
GIt=2.7216*10^8 # Ncm2      # väändejäikus
baasi0=1.0*10^10   # baasjäikus
EIw=7.5600*10^12 #Ncm4 # kooldumisjäikus
#######
spvF=yspWGvfmhvI(baasi0,l,GIt,EIw);  # hõreda laiendatud ülekandemaatriksi arvutus
vB=yzWGMx(baasi0,l,l,a,Mx,GIt,EIw); # koormusvektori arvutus
IIv=1;
IJv=1;
vB=vB.*(-1);
spA=spInsertBtoA(spA,IIv,IJv,spvF);     # astab ülekandemaatriksi võrrandisüsteemi
Bvb=InsertBtoA(Bvb,NNK,1,IIv,1,vB,4,1); # astab koormusvektori vabaliikmesse
####### Rajatingimused 
k2=GIt/EIw;
spA=spSisestaArv(spA,5,1,1);  # $theta_A$ - väändenurk
spA=spSisestaArv(spA,6,3,1);  # $B_{L}$ -bimoment 
spA=spSisestaArv(spA,7,5,1);  # $theta_A$ - väändenurk
spA=spSisestaArv(spA,8,7,1);  # $B_{L}$ -bimoment 
# Bvb(8,1)=10000.0 ;
#  vastavd  vabaliikmed Bvb on juba nullistatud
#
spA_rank = sprank(spA) # võrrandisüsteemi astak
A=full(spA)
Bvb
# Teine märgikokkulepe.Väändnurk skaleeritud
Z=spA\Bvb; # lahendab võrrandisõteemist spA*X=Bvb
Z

#
## OK lk 52 kõik klapib
  
NEARV=1;
#baasi0=1.0
#for i=1:NEARV
#%siireVardaA(i,1:2)
AlgPar(1:4,1)=[Z(1,1)/baasi0  Z(2,1) Z(3,1)  Z(4,1)];
#endfor
#      
#AlgParm=AlgPar;  
disp('skaleerimata väändenurk skaleerimata väändenurk varda lõpus ')
theta=Z(5,1)/baasi0 # skaleerimata väändenurk varda lõpus
spA  # hõre laiendatud ¨ulekandemaatriks
disp(' Algparameetrid - AP ') 
   AP=AlgPar(:,1)
baasi0=1.0
Nmitmeks=4
xx=0;
xsamm=l/Nmitmeks;
for ij=1:Nmitmeks+1  # 5 - displacements and forces at x=0.0   
Xloikes(ij,1)=xx;
   vvF=ylWGfhlin(baasi0,l,xx,GIt,EIw);
# vvF=ylTfhlin(1.0,xx,GAr,EI);               # ylfhlin(1.0,xx,EA,GAr,EI);  
  vvB=yzWGMx(baasi0,l,xx,a,Mx,GIt,EIw); 
#  vvB=ESTtalaKrmus(1.0,xx,Li,Fjoud,qkoormus,EI);
  
 Fvv(1:4,ij)=vvF*AP+vvB;
 Fvv(5,ij)=Fvv(2,ij)+Fvv(4,ij);
xx=xx+xsamm;
endfor 
#Fvv(1,:)
Fvv;
kt=sqrt(GIt/EIw);
k=kt
  nu=l*sqrt(GIt/EIw);
 nul=nu/l;
 lnu=l/nu;
# nux=nul*x;
  nux=nul*l;
 ch=cosh(nux);
 sh=sinh(nux);
#
  
  suurused=[' theta -'; 
            '    Tt -'; 
            '     B -'; 
            '    Tw -'; 
            '  Tsum -']; 
#disp(' size(suurused) ');            
  size(suurused);          

#for i=1:5
disp(sprintf('%9s %9.2f    %9.2f     %9.2f    %9.2f    %9.2f', ' x= ',  Xloikes(1,1), Xloikes(2,1), Xloikes(3,1), Xloikes(4,1),Xloikes(5,1)))
#endfor
 # 
#
for i=1:5 
disp(sprintf('%8s %12.3e %12.3e %12.3e %12.3e %12.3e',suurused(i,:), Fvv(i,1), Fvv(i,2), Fvv(i,3), Fvv(i,4), Fvv(i,5)))
endfor
 # 
disp('  ')  
X=l/2-eps1;
SFv=SisejoudWGMxpunktis(l,X,a,AlgPar,Mx,GIt,EIw,suurused); 
 # 
disp('  ') 
#%%%%%%%%%%%%%%%%%%%%%%%%%%%
#figure(1)
#spy(spA)
#title('spy(spA) - the sparse matrix spA(8,8) non zero elements [31%] ')
#refresh  
##print('-dfig','-landscape','spESTNaide2_2a_sparse_matrix.fig');  
refresh
disp('-----------------------------------')
disp('  ')
disp(' The END ')
#%%%%%%%%%%%%%%%  
 #  
 MiMaSg=minMaxGrfVnurk(Fvv);
disp('  ') 
Nmitmeks=800
xx=0;
xsamm=l/Nmitmeks;
for ij=1:Nmitmeks+1  # 5 - displacements and forces at x=0.0   
Xloikes(ij,1)=xx;
   vvF=ylWGfhlin(baasi0,l,xx,GIt,EIw);
# vvF=ylTfhlin(1.0,xx,GAr,EI);               # ylfhlin(1.0,xx,EA,GAr,EI);  
  vvB=yzWGMx(baasi0,l,xx,a,Mx,GIt,EIw);
  
 Fvv(1:4,ij)=vvF*AP+vvB;
 Fvv(5,ij)=Fvv(2,ij)+Fvv(4,ij);
 Fvv(6,ij)=xx;
xx=xx+xsamm;
ViimaneX=ij;
endfor 
# 
# %%%%%%%%%%%%%%%%%%%%%%%%%%%
ViimaneX=ij
disp('  ')  
X=l/2-eps2;
SFv=SisejoudWGMxpunktisGr(l,X,a,AlgPar,Mx,GIt,EIw,suurused); 
 # 
disp('  ') 
ViimaneX=ViimaneX+1;
Fvv(1:6,ViimaneX)=SFv(1:6,1);
FvvSize=size(Fvv);
FvvTus=Fvv';
FvvT=sortrows(FvvTus,6);
FvvT(:,6);
Fvv=FvvT';
#%%%%%%%%%%%%%%%%%%%%%%%%%%%
figure(1)
spy(spA)
title('spy(spA) - the sparse matrix spA(8,8), 20 non zero elements [31%] ')
refresh  
#print('-dfig','-landscape','spESTNaide2_1_sparse_matrix.fig');  
refresh
disp('-----------------------------------')
disp('  ')
#disp(' The END ')
#%%%%%%%%%%%%%%% 
 % 
disp('  ') 
disp('  ') 
%%%%%%%%%%%
Fvv5(1:5,:)=Fvv(1:5,:);
  MiMaTa=minMaxGrfVnurk(Fvv5);
 %   MiMaTa=minMaxGrfVnurk(Fvv);
  	maxthetagr=MiMaTa(1,1);
	minthetagr=MiMaTa(1,2);
	maxTtgr=MiMaTa(2,1);
	minTtgr=MiMaTa(2,2);
	maxBwgr=MiMaTa(3,1);
	minBwgr=MiMaTa(3,2);
	maxTwgr=MiMaTa(4,1);
	minTwgr=MiMaTa(4,2);        
	maxTsumgr=MiMaTa(5,1);
	minTsumgr=MiMaTa(5,2); 
%%%%%%%%%
L=l;
JA=-0.1*L;
JL=l+0.1*L;
jaotT=L/8;
XL1=L;
Xg=[0.0; XL1];
Yg=[0.0; 0.0];
%%%%%%%%%
XloikesT=FvvT(:,6)';
 FvvT=Fvv; 
%%
 XA0=XloikesT(1,1);
 YA0=0.0;

FvvSize=size(Fvv);
ViimaneX=FvvSize(1,2);
%%%%%%%%% 
 YAY=FvvT(1,1);
 YLY=FvvT(1,ViimaneX);
 X0=[0.0; 0.0];
 Y0=[0.0;  YAY];
 YAY=FvvT(1,1);
 XL=[L; L];
 YL=[0.0; YLY];
 
%%Xgraf= size(XloikesT)
%%Ygraf=size(FvvT)

%%%%%%%%%   
figure(2)
%
 plot( X0,Y0,"3",XloikesT(1,:),FvvT(1,:),"3",Xg,Yg,"3",XL,YL,"3")
grid("on")
% http://www.network-theory.co.uk/docs/octave3/octave_160.html
%grid("off")
%minthetagr-0.00002
ax = gca ();
set (ax, "outerposition", [ 0  0  1 0.55])
%minthetagr-150maxthetagr+150
 axis([JA JL  minthetagr-0.001 maxthetagr+0.002], "ij")  % telgede piirid [Xalg Xlõpp Zalg Zlõpp]
 ax=gca();
   xticks(1,:)=0.0:jaotT:L;
 set (ax, "xtick", xticks)
title('V66ndenurk theta ',"fontsize",10)
        xlabel('x [m]',"fontsize",10)
	ylabel('theta [rad]',"fontsize",10)
%
%%print('-deps','V66ndenurk_theta.eps');
%%print('-dpng','V66ndenurk_theta.png');
%print('-dfig','-landscape','-mono','-solid','V66ndenurk_theta.fig');
%%print('-dfig','-solid','-landscape','-color','V66ndenurk_theta.fig');
%%http://www.gnu.org/software/octave/doc/interpreter/Printing-Plots.html
refresh
%title('Kooldumisv66ndemoment Tw ',"fontsize",10)
%Vabav66ndemoment
%Koguv66ndemoment
%%%%%%%%%    
%%%%%%%%% 
 YAY=FvvT(2,1);
 YLY=FvvT(2,ViimaneX);
 X0=[0.0; 0.0];
 Y0=[0.0;  YAY];
 YAY=FvvT(2,1);
 XL=[L; L];
 YL=[0.0; YLY];
%%%%%%%%%  
%figure(3)  
 figure(3)
%
 plot( X0,Y0,"3",XloikesT(1,:),FvvT(2,:),"3",Xg,Yg,"3",XL,YL,"3")
grid("on")
% http://www.network-theory.co.uk/docs/octave3/octave_160.html
%grid("off")
%
ax = gca ();
set (ax, "outerposition", [ 0  0  1 0.55])
%
 axis([JA JL minTtgr-2000 maxTtgr+2000], "ij")  % telgede piirid [Xalg Xlõpp Zalg Zlõpp]
 ax=gca();
   xticks(1,:)=0.0:jaotT:L;
 set (ax, "xtick", xticks)
title('Vabav66ndemoment Tt ',"fontsize",10)
        xlabel('x [m]',"fontsize",10)
	ylabel('Tt [kNcm]',"fontsize",10)
%
%%print('-deps','Vabav66ndemoment_Tt.eps');
%%print('-dpng','Vabav66ndemoment_Tt.png');
%print('-dfig','-landscape','-mono','-solid','Vabav66ndemoment_Tt.fig');
%%print('-dfig','-solid','-landscape','-color','Vabav66ndemoment_Tt.fig');
%%http://www.gnu.org/software/octave/doc/interpreter/Printing-Plots.html
refresh
%title('Kooldumisv66ndemoment Tw ',"fontsize",10)
%Vabav66ndemoment

%%%%%%%%% 
 YAY=FvvT(3,1);
 YLY=FvvT(3,ViimaneX);
 X0=[0.0; 0.0];
 Y0=[0.0;  YAY];
 XL=[L; L];
 YL=[0.0; YLY];
%%%%%%%%% 
figure(4)
%
 plot( X0,Y0,"3",XloikesT(1,:),FvvT(3,:),"3",Xg,Yg,"3",XL,YL,"3")
grid("on")
% http://www.network-theory.co.uk/docs/octave3/octave_160.html
%grid("off")
%
ax = gca ();
%set (ax, "outerposition", [ 0  0  1 0.4])
set (ax, "outerposition", [ 0  0  1 0.55])
%
 axis([JA JL minBwgr-300000 maxBwgr+300000], "ij")  % telgede piirid [Xalg Xlõpp Zalg Zlõpp]
 ax=gca();
   xticks(1,:)=0.0:jaotT:L;
 set (ax, "xtick", xticks)
title('Bimoment Bw ',"fontsize",10)
        xlabel('x [m]',"fontsize",10)
	ylabel('B [kNcm2]',"fontsize",10)
%
%%print('-deps','Bimoment_Bw.eps');
%%print('-dpng','Bimoment_Bw.png');
%print('-dfig','-landscape','-mono','-solid','Bimoment_Bw.fig');
%%print('-dfig','-solid','-landscape','-color','Bimoment_Bw.fig');
%%http://www.gnu.org/software/octave/doc/interpreter/Printing-Plots.html
refresh
% 
%%%%%%%%% 
 YAY=FvvT(4,1);
 YLY=FvvT(4,ViimaneX);
 X0=[0.0; 0.0];
 Y0=[0.0;  YAY];
 YAY=FvvT(4,1);
 XL=[L; L];
 YL=[0.0; YLY];
%%%%%%%%%  
figure(5)
%
 plot( X0,Y0,"3",XloikesT(1,:),FvvT(4,:),"3",Xg,Yg,"3",XL,YL,"3")
grid("on")
% http://www.network-theory.co.uk/docs/octave3/octave_160.html
%grid("off")
%
ax = gca ();
set (ax, "outerposition", [ 0  0  1 0.55])
%
 axis([JA JL minTwgr-5000 maxTwgr+5000], "ij")  % telgede piirid [Xalg Xlõpp Zalg Zlõpp]
 ax=gca();
   xticks(1,:)=0.0:jaotT:L;
 set (ax, "xtick", xticks)
title('Kooldumisv66ndemoment Tw ',"fontsize",10)
        xlabel('x [m]',"fontsize",10)
	ylabel('Tw [kNcm]',"fontsize",10)
%
%%print('-deps','Kooldumisv66ndemoment_Tw.eps');
%%print('-dpng','Kooldumisv66ndemoment_Tw.png');
%print('-dfig','-landscape','-mono','-solid','Kooldumisv66ndemoment_Tw.fig');
%%print('-dfig','-solid','-landscape','-color','Kooldumisv66ndemoment_Tw.fig');
%%http://www.gnu.org/software/octave/doc/interpreter/Printing-Plots.html
refresh
%title('Kooldumisv66ndemoment Tw ',"fontsize",10)
%Vabav66ndemoment

%%%%%%%%% 
 YAY=FvvT(5,1);
 YLY=FvvT(5,ViimaneX);
 X0=[0.0; 0.0];
 Y0=[0.0;  YAY];
 YAY=FvvT(5,1);
 XL=[L; L];
 YL=[0.0; YLY];
%%%%%%%%%   
figure(6)
%
 plot( X0,Y0,"3",XloikesT(1,:),FvvT(5,:),"3",Xg,Yg,"3",XL,YL,"3")
grid("on")
% http://www.network-theory.co.uk/docs/octave3/octave_160.html
%grid("off")
%
ax = gca ();
set (ax, "outerposition", [ 0  0  1 0.55])
%
 axis([JA JL minTwgr-6000.0 maxTwgr+6000.0], "ij")  % telgede piirid [Xalg Xlõpp Zalg Zlõpp]
 ax=gca();
   xticks(1,:)=0.0:jaotT:L;
 set (ax, "xtick", xticks)
title('Koguv66ndemoment Tsum ',"fontsize",10)
        xlabel('x [m]',"fontsize",10)
	ylabel('Tsum [kNcm]',"fontsize",10)
%
%%print('-deps','Koguv66ndemoment_Tsum.eps');
%%print('-dpng','Koguv66ndemoment_Tsum.png');
%print('-dfig','-landscape','-mono','-solid','Koguv66ndemoment_Tsum.fig');
%%print('-dfig','-solid','-landscape','-color','Koguv66ndemoment_Tsum.fig');
%%http://www.gnu.org/software/octave/doc/interpreter/Printing-Plots.html
refresh
%title('Kooldumisv66ndemoment Tw ',"fontsize",10)
%Vabav66ndemoment
%Koguv66ndemoment
%%%%%%%% 
  
  
  
  
  
  
  
   
   
   
   
  
   
   
   
   








