##  LaheBeamDFI.m
function AlgParm=LaheBeamDFI(baasi0,Ntoerkts,esQkoormus,esFjoud,sSolmF,tsolm,tSiire,krdn,selem)
##
## 
##disp('==================================================================')
##disp(' LaheBeamDFI.m     ')
##disp(' The compiler of a basic-, compatibility-, joint equilibrium-, ')
##disp(' side conditions- and support restrtictions sparse equations.  ')
##disp(' Solving system of sparse equations. Find the initial parmeter vectors ')
##disp(' for elements displacements, forces and the support reactions.  ')
##disp(' OUTPUT: AlgParm -- the initial parmeter vector for elements. ')
## 
## GNU octave version >= 3.0.x
##=========================================================================
## This Program is written by Andres Lahe,   2010-07-05
##                    e-mail: andres.lahe@ttu.ee
## LAST MODIFIED: Andres Lahe,   2013-07-08
## Copyright (c)  2010 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
##=========================================================================
##
## baasi0 - scaling multiplier for the displacements,
## Ntoerkt - the number of the support reactions,
## esQkoormus - the uniformly distributed load in local coordinate z direction
##              esQkoormu(LoadsqONelement,3,ElementideArv),
## esFjoud - the point load in local coordinate z direction
##           esFjoud(LoadsF_on_Element,2,ElementideArv),
## sSolmF - the node forces in global  coordinates
##          sSolmF(2,1,SolmedeArv),
## tsolm - the restrtictions on the support displacements 
## tSiire - the support shift 
##          tSiire(2,1,SolmedeArv)
## krdn - the nodal coordinates
## selem - the topology of elements 
##         selem(ElementideArv,1:16)
##
#if nargin != 9
if ~(nargin==9)
error(' function LaheBeamDFI have wrong number of input arguments!') 
end
#
 global spA
 global X
 global siireVardaA
#disp('---  ') 
#%selem=[selemjl(:,1:23)];
EARV=size(selem);
NEARV=EARV(1,1);
SARV=size(krdn);
NSARV=SARV(1,1);
NNK=8*NEARV+Ntoerkts;
#==========
TSARV=size(tsolm);
NTSARV=TSARV(1,1);
##%%%%%%%%%%%
ToeSidemeteArv=0;
for i=1:NTSARV
sJrN(i)=i;
  for j=2:3
     if (tsolm(i,j) == 1)
       ToeSidemeteArv=ToeSidemeteArv+1;
     endif 
  endfor
endfor
#
#%%%%%%%%%%%%%%%%%%%%%%
selemjl=selem;
#1
lvarras=VardaPikkusT(NSARV,NEARV,krdn,selem);
#lvarras
disp('   ')
NNK=8*NEARV+Ntoerkts;
#%%%%%%%%%%%
#
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#
#  Nr;  daf''s at the end:  w, fi, Q, M; Node number; shear-, moment hinge.
#disp(' Nr;  daf''s at the end:  w, fi,, Q, M; Node number; shear-, moment hinge ')
#
for i=1:NEARV
selemlopp(i,1)=i;
selemlopp(i,2:5)=selemjl(i,1:4);
selemlopp(i,6)=[selemjl(i,11)];
selemlopp(i,7:8)=[selemjl(i,13:14)];
endfor
#
#selemlopp
#%%%%%%%%%%%
# Nr;  daf''s at the beginning:  w, fi, Q, M; Node number; shear-, moment hinge.
#disp(' Nr;  daf''s at the beginning:  w, fi, Q, M; Node number; shear-, moment hinge ')
#
for i=1:NEARV
selemaluses(i,1)=i;
selemaluses(i,2:5)=selemjl(i,5:8);
selemaluses(i,6)=[selemjl(i,12)];
selemaluses(i,7:8)=[selemjl(i,15:16)];
endfor
#
%selemaluses
#%%%
%elemvabsolm
for i=1:NEARV
elemvabsolm(i,1:8)=selemlopp(i,1:8);
endfor
#
for i=1:NEARV
elemvabsolm(NEARV+i,1:8)=selemaluses(i,1:8);
endfor
#%%%
elemvabsolm;
AB=sortrows(elemvabsolm(:,6));
#AB;
#
##[~,idx] = sortrows(elemvabsolm(:,6));
[tmp,idx] = sortrows(elemvabsolm(:,6));
ABB=elemvabsolm(idx,:);
disp('============================================================= ')
disp('Element number; DaF numbers; Node number; shear, moment hinge')
disp('------------------------------------------------------------- ')
disp('Element number  w fi  Q M   Node  Q    M  hinge -  true=1     ')
disp('------------------------------------------------------------- ')
ABB=ABB
#--------- 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('  B=zeros(NNK,1);  ')
B=zeros(NNK,1);   % nullistan võrrandisüsteemi vabaliikmed
disp('   ')
%%%%%%%%%%%%%%%%%%%%%%%%%

#------------ The basic equations of beam ----------------
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
disp('----- Writing basic equations of beam with the transfer matrix ----  ')
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
#  http://digi.lib.ttu.ee/opik_eme/Ehitusmehaanika.pdf#page=397   page=693
##%%%%%%%%%%%%%%%%%%%%%%%
IIv=0;
IJv=0;
#
for i=1:NEARV
qkoormus=zeros(1,3);
Fjoud=zeros(1,2);
spvF=zeros(4,8);
vB=zeros(4,1);
vFz=zeros(4,1);
#Sj=zeros(4,1);
krda=i;
EI=selem(i,9);
GAr=selem(i,10);
Li=lvarras(i,1);
#%qx=qxZ(i,1); # qxZ(2,1)
#AL#qz=qzZ(i,1);
#AL#aLx=aLXx(i,1);
#AL#Fz=FZz(i,1);
#
Fjoud=esFjoud(:,1:2,i);
qkoormus=esQkoormus(:,1:3,i);
##%Fx=FZx(i,1);
spvF= yspTlvfmhvI(baasi0,Li,Li,GAr,EI); 
#vB=yzThqz(baasi0,Li,qz,EI);  
#vFz=yzTfzv(baasi0,Li,aLx,Fz,EI); 
#vB=vB+vFz;
##vB=ESTtalaKrmusVr(baasi0,Li,Fjoud,qkoormus,EI);
##vB=ESTtalaKrmus(baasi0,Li,Li,Fjoud,qkoormus,EA,EI);
vB=ESTtalaKrmus(baasi0,Li,Li,Fjoud,qkoormus,EI);
#disp('real 183  ') 
IIv=krda*4-3;
IJv=krda*8-7;
spA=spInsertBtoA(spA,IIv,IJv,spvF);
B=InsertBtoA(B,NNK,1,IIv,1,vB,4,1);
#
endfor
#
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#
vorrandeidB=size(spA);
rows_of_basic_equations=vorrandeidB(1,1);
cols_of_basic_equations=vorrandeidB(1,2);
non_zero_elements_in_the_basic_equations=nnz(spA);

spRidu=vorrandeidB(1,1);
pohivorrandiXe=vorrandeidB(1,2);  # row 651
disp('  ')
disp('--------The basic equations are inserted into spA  --------  ')
rows=sprintf('rows_of_basic_equations: %d', spRidu)
col=sprintf('cols_of_basic_equations: %d', pohivorrandiXe)
spA_nnz=sprintf('non_zero_elements_in_the_basic_equations: %d', nnz(spA))
nnzB=spA_nnz;
disp('  ')
Nr1=spRidu+1;
#
#disp('  ')
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%spA
%B
#%%%%%%%%%%%
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  VarrasteArv=NEARV;
  SolmedeArv=NSARV;
#%%%%
    Vord='=';
    Kom=',' ;  
    Semi=';';
   #Tyhk='  ';     
 Arv1=num2str(1);
 Arv1m=num2str(-1); 
 Arv2=num2str(2);
% Arv3=num2str(3); 

%disp(' rida 227 ')
%spA
%spA=sparse(NNK,NNK)
%pause
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
disp(' Compatibility equations of displacements at nodes ')
From_row=sprintf('Compatibility equations begin from row: %d', Nr1)
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
#
for i=1:NSARV
#disp('---------------------------------------------------------')
Solm=i;
Node=Solm
%%%
%AB
%ABB
AVS=VardadSolmesT(NSARV,NEARV,Solm,AB,ABB);
#
#disp('---------------------------------------------------------') 
#
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
NEARV2=NEARV*2;
#
SpTM2x2xz=sparse(2,2);
SpTM1x0xz=sparse(1,2);
SpTM0x1xz=sparse(1,2);
#
AVSdim=size(AVS) ;

if AVSdim(1,1) > 1

    for j=2:AVSdim(1,1)
    
       Varras1=AVS(1,1);    
       VarrasN=AVS(j,1);
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#
   str1=num2str(Nr1) ;
%%   AVS
% AVS18=AVS(1,8)
%   AVSJ8=AVS(j,8)
   AVSsum=AVS(1,8)+ AVS(j,8);
#
    if AVSsum == 0
      SpTM2x2xz=SpTeisendusUhikMaatriks2x2(1);
      TAVS2=eval(sprintf('%3i',AVS(1,2)));  
      TVS2=int2str(TAVS2);
      str2=TVS2;    
      str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM2x2xz',')',Semi];
      cmd = sprintf(str)
      eval(cmd)  
      SpTM2x2xz=SpTeisendusUhikMaatriks2x2(1);
      SpTM2x2xzm=-SpTM2x2xz;
# disp(' real 280 -SpTM2x2xz ')     
      TAVS2=eval(sprintf('%3i',AVS(j,2)));  
      TVS2=int2str(TAVS2);
      str2=TVS2;
      str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM2x2xzm',')',Semi];
      cmd = sprintf(str)
      eval(cmd);
      Nr1=Nr1+2; 
    else    
# 
      SpTM1x0xz=SpTeisendusUhikMaatriks1x0(1);
#      SpTM2x2;
      str1=num2str(Nr1);    
      TAVS2=eval(sprintf('%3i',AVS(1,2)));  
      TVS2=int2str(TAVS2);
      str2=TVS2;    
      str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM1x0xz',')',Semi];
      cmd = sprintf(str)
      eval(cmd); 
      SpTM1x0xz=SpTeisendusUhikMaatriks1x0(1);
      SpTM1x0xzm=-SpTM1x0xz; 
      TAVS2=eval(sprintf('%3i',AVS(j,2)));  
      TVS2=int2str(TAVS2);
      str2=TVS2;
      str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM1x0xzm',')',Semi];
      cmd = sprintf(str)
      eval(cmd);        
      Nr1=Nr1+1; 
    endif
# 
  endfor  #% for j=2:AVSdim(1,1)                   
 endif #%if AVSdim(1,1) > 1                                       
endfor  #% for i=1:NSARV
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
disp('----- ')
vorrandeidD=size(spA);
#
disp(' ')
spA_rows=vorrandeidD(1,1)
spA_rowsD=spA_rows;
spA_cols=vorrandeidD(1,2)
#
nnzD=nnz(spA);
spA_nnz=sprintf('non_zero_elements_in_spA: %d', nnz(spA))
disp('  ')
disp('-----Compatibility equations of displacements are inserted into spA  ----  ')
disp('  ')
compatibility_equations_rows=spA_rowsD-rows_of_basic_equations
non_zero_elements_in_compatibility_equations= nnzD-non_zero_elements_in_the_basic_equations
#nnzB=nnz(spA)
disp('  ')
#%
%disp(' rida 332 ')
%spA=sparse(NNK,NNK)
%pause
#%pohivorrandiXe
toemuutuja=pohivorrandiXe+1; # row 336  # row 651
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
disp(' Joint equilibrium equations at nodes ')
From_rows=sprintf('Joint equilibrium equations begin from row: %d', Nr1)
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
#
Nr1;
koolonV=':';
#
tsolmLaiend=zeros(NSARV,3); # Nodes with support restrtictions
  for i=1:NSARV
      tsolmLaiend(i,1)=i; 
    for j=1:NTSARV;
     if (tsolmLaiend(i,1)-tsolm(j,1)) == 0
      tsolmLaiend(i,2:3)=tsolm(j,2:3);            %% Siin viga ?
     endif 
    endfor
 endfor
#
##tsolm
 tsolmLaiend;
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:NSARV  % sõlmede arv
#     uvfiLaiend(i)=tsolmLaiend(i,2)+tsolmLaiend(i,3);
#  disp('---------------------------------------------------------')
  Solm=i;
  Node=Solm;
  str1=num2str(Nr1); 
  strSolm=num2str(Solm);
  AVS=VardadSolmesT(NSARV,NEARV,Solm,AB,ABB)
      ##  AVS 
  wfiLaiend(Node,1)=tsolmLaiend(Node,2)+tsolmLaiend(Node,3);  #  == support restrtictions
#
#  disp('---------------------------------------------------------') 
  vorrandeid=size(spA);
  spRidu=vorrandeid(1,1); 
%% 
#%%%%%%%%%%%%%%%%%%
  NEARV2=NEARV*2;
  AVSdim=size(AVS) ;
  ElArvSolmes=AVSdim(1,1);
#%%%%%%%%%%%%%%%%%%%  
## Nodes with support restrtiction. wfiLaiend(Node,1) - the number of restrtiction.  
## 
   if (wfiLaiend(Node,1) != 0)
       # with restrtictions
      if (ElArvSolmes -2 ) == 0 
 #         disp(' 2 elements ')
 #         Hinge?  AVS(1,8)  0 või 1
 #         AVS(1 ,1)  ---  AVS(1 ,4) - jouNr
 #         AVS(2 ,1)   --- AVS(2 ,4) - jouNr
            if (AVS(1,8) -1) == 0 
#                   disp(' Hinge yes 389 ')
#                   Node=Node
                for j=1:AVSdim(1,1)   
                  VarrasN=AVS(j,1);
                  str1=num2str(Nr1) ;
                  AVStas= AVS(j,8); 
                  SpTM1x0xz=SpTeisendusUhikMaatriks1x0(1);
                  TAVS2=eval(sprintf('%3i',AVS(j,4)));  
                  TVS2=int2str(TAVS2);
                  str2=TVS2;    
                  str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM1x0xz',')',Semi];
                  cmd = sprintf(str)
                  eval(cmd) 
#         disp('      Nodal forces at the node')
                  strSolmB=num2str(Solm);
                  strB1=str1;
                  NrB2=Nr1+2;
                  strB2=num2str(NrB2);
                  Arv123=Arv2;
                endfor 
#                 disp('      Nodal forces at the node')
                  strB=['B(',strB1, koolonV, strB2, Kom,Arv1,')','=','sSolmF(', Arv1, koolonV, Arv2, Kom, strSolmB,')',Semi];
                  cmd = sprintf(strB)
                  eval(cmd) 
                  Nr1=Nr1+1;
             else 
#                    disp(' no hinges 402')
                    Node=Node
                for j=1:AVSdim(1,1)   
                  VarrasN=AVS(j,1);
                  str1=num2str(Nr1) ;
                  AVStas= AVS(j,8);   
                  SpTM2x2=SpTeisendusUhikMaatriks2x2(1);;  # rida 420
                  TAVS2=eval(sprintf('%3i',AVS(j,4)));  
                  TVS2=int2str(TAVS2);
                  str2=TVS2;    
                  str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM2x2',')',Semi];
                  cmd = sprintf(str)
                  eval(cmd) 
 #        disp('  770    Nodal forces at the node')
                  strSolmB=num2str(Solm);
                  strB1=str1;
                  NrB2=Nr1+1;
                  strB2=num2str(NrB2);
                  Arv123=Arv2;
                endfor # for j=1:AVSdim(1,1) 
#                 disp('      Nodal forces at the node')
                  strB=['B(',strB1, koolonV, strB2, Kom,Arv1,')','=','sSolmF(', Arv1, koolonV, Arv2, Kom, strSolmB,')',Semi];
                  cmd = sprintf(strB)
                  eval(cmd) 
                  Nr1=Nr1+2; 
             endif #  if  (AVS(1,8) -1) == 0 
       else
%           disp(' 1 element  785')
%           Kas liigend?  AVS( ,8)  0 või 1
              if (AVS(1,8) -1) == 0 
                   disp(' Hinge yes 792')
                   
            else 
#                   disp(' No hinges  235')
                   Node=Node
              for j=1:AVSdim(1,1)   
                   VarrasN=AVS(j,1);
                   str1=num2str(Nr1) ;
                   AVStas= AVS(j,8); 
#               disp('  SpTM2x2 rida 453 ')   
                   SpTM2x2=SpTeisendusUhikMaatriks2x2(1);;  # rida 454 SpTM2x2=SSpTM2x2xz=SpTei
                   TAVS2=eval(sprintf('%3i',AVS(j,4)));  
                   TVS2=int2str(TAVS2);
                   str2=TVS2;    
                   str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM2x2',')',Semi];
                   cmd = sprintf(str)
                   eval(cmd)  
#         disp('      Nodal forces at the node')
                   strSolmB=num2str(Solm);
                   strB1=str1;
                   NrB2=Nr1+1;
                   strB2=num2str(NrB2);
                   Arv123=Arv2;
                endfor # for j=1:AVSdim(1,1)
 #                 disp('      Nodal forces at the node')
                   strB=['B(',strB1, koolonV, strB2, Kom,Arv1,')','=','sSolmF(', Arv1, koolonV, Arv2, Kom, strSolmB,')',Semi];
                   cmd = sprintf(strB)
                   eval(cmd) 
                   Nr1=Nr1+2;  
             endif # if (AVS(1,8) -1) == 0 
       endif  # if (ElArvSolmes -2 ) == 0
       
 # substitute supports into the equations 
       
          if (wfiLaiend(Node,1) -2 ) == 0 
#               Node
               VarrasN=AVS(1,1);   # ???????????
##             disp('Q M restrtictions ')
               SpTM2x2xz=SpTeisendusUhikMaatriks2x2(1);
               SpTM2x2xzm=-SpTM2x2xz;
#   disp(' real 484  esimesel toel toereaktsioonid paika ')         
#   disp(' real 485 -SpTM2x2xz ')             
               toemuutuja;
               TAVS2=toemuutuja;
               TVS2=int2str(TAVS2);
               str2=TVS2;    
               str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM2x2xzm',')',Semi];
               cmd = sprintf(str) 
               eval(cmd)
               toemuutuja=toemuutuja+2;
             
            else
               if (tsolmLaiend(Node,2) == 1) 
                  Node;
#                disp('Q restrtictions ') 
                  % VarrasN=AVS(1,1)   # ???????????
#                disp('real 500  -SpTM1x0xz  on tugi 3 , 6, 7 Z suunaline')                
                  SpTM1x0xz=SpTeisendusUhikMaatriks1x0(1);
                  SpTM1x0xzm=-SpTM1x0xz;
#                   str1
#                  toemuutuja
                  TAVS2=toemuutuja; 
                  TVS2=int2str(TAVS2);
                  str2=TVS2;    
                  str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM1x0xzm',')',Semi];
                  cmd = sprintf(str) 
                  eval(cmd)
                  toemuutuja=toemuutuja+1;                                    
              endif
           
               if (tsolmLaiend(Node,3) == 1) 
                  Node; 
#                 VarrasN=AVS(1,1);   # ??????????
                  disp('M restrtictions ')
                  SpTM0x1xz=SpTeisendusUhikMaatriks0x1(1);
                  SpTM0x1xzm=-SpTM0x1xz;
                  TAVS2=toemuutuja;
                  TVS2=int2str(TAVS2);
                  str2=TVS2;    
                  str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM0x1xzm',')',Semi];
                  cmd = sprintf(str) 
                  eval(cmd)
                  Nr1=Nr1+1;
                  toemuutuja=toemuutuja+1;
                endif   
             endif ##(wfiLaiend(Node,1) -2 ) == 0
     else  
        #  
#        disp(' No suppors. No restrtictions ')
              Node=Node
        if (ElArvSolmes -2 ) == 0 
        
 #         disp(' 2 elements ')
 #         Hinge?  AVS(1,8)  0 või 1
 #         AVS(1 ,1)  ---  AVS(1 ,4) - jouNr
 #         AVS(2 ,1)   --- AVS(2 ,4) - jouNr
 
              if (AVS(1,8) -1) == 0 
#                   disp(' hinge yes 543 ')
#                   Node=Node
                for j=1:AVSdim(1,1)   
                   VarrasN=AVS(j,1);
#                    Node=Node
                   str1=num2str(Nr1) ;
                   AVStas= AVS(j,8);   
##               %%   SpTM2x2=;  # row 535 SpTM2x2xz=SpTeisendusUhikMaatriks2x2(1);
                   SpTM1x0xz=SpTeisendusUhikMaatriks1x0(1);
                   TAVS2=eval(sprintf('%3i',AVS(j,4)));  
                    TVS2=int2str(TAVS2);
                    str2=TVS2;    
                    str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM1x0xz',')',Semi];
                    cmd = sprintf(str)
                    eval(cmd)  
#         disp('      Nodal forces at the node')
                    strSolmB=num2str(Solm);
                    strB1=str1;
                    NrB2=Nr1;    #  +1;         ####### Siin on  NrB2=Nr1+2; peab olema  NrB2=Nr1+1;
#   disp(' real 561     strB2=num2str(NrB2) ')             
                   strB2=num2str(NrB2);  ## ??????????????
                   Arv123=Arv2;
                endfor # for j=1:AVSdim(1,1)   
#                 disp('      Nodal forces at the node  938  ')
                    strB=['B(',strB1, koolonV, strB2, Kom,Arv1,')','=','sSolmF(', Arv1, koolonV, Arv1, Kom, strSolmB,')',Semi];
                %%    strB=['B(',strB1, koolonV, strB2, Kom,Arv1,')','=','sSolmF(', Arv1, koolonV, Arv2, Kom, strSolmB,')',Semi];
                    cmd = sprintf(strB)
                    eval(cmd) 
                    Nr1=Nr1+1;
              else 
#                    disp(' no hinges 572')
#                    disp(' hinge yes 573 ')
                  Node=Node
                for j=1:AVSdim(1,1)
                      Node=Node
                    VarrasN=AVS(j,1);
                    str1=num2str(Nr1);
                    AVStas= AVS(j,8);   
                    SpTM2x2=SpTeisendusUhikMaatriks2x2(1);   # rida 580
                    TAVS2=eval(sprintf('%3i',AVS(j,4)));  
                    TVS2=int2str(TAVS2);
                    str2=TVS2;    
                    str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM2x2',')',Semi];
                    cmd = sprintf(str)
                    eval(cmd)  
 #        disp('  770    Nodal forces at the node')
                    strSolmB=num2str(Solm);
                    strB1=str1;
                    NrB2=Nr1+1;
                    strB2=num2str(NrB2);
                    Arv123=Arv2;
                 endfor # for j=1:AVSdim(1,1)
#              disp('      Nodal forces at the node')
                   strB=['B(',strB1, koolonV, strB2, Kom,Arv1,')','=','sSolmF(', Arv1, koolonV, Arv2, Kom, strSolmB,')',Semi];
                   cmd = sprintf(strB)
                   eval(cmd) 
                   Nr1=Nr1+2; 
               endif # if (AVS(1,8) -1) == 0
       else
%           disp(' 1 elemen  785')
%           Hinge?  AVS( ,8)  0 või 1
                if (AVS(1,8) -1) == 0 
                   disp(' hinge yes 792')
                else 
#                    disp(' no hinge  594')
         for j=1:AVSdim(1,1)   
                  VarrasN=AVS(j,1);
                  str1=num2str(Nr1);
                  AVStas= AVS(j,8);   
#             disp(' rida 611   pTM2x2=SSpTM2x2xz= ')  
                  SpTM2x2=SpTeisendusUhikMaatriks2x2(1); ## SpTM2x2=SSpTM2x2xz=SpTeisendusUhikMaatriks2x2
                 
                  TAVS2=eval(sprintf('%3i',AVS(j,4)));  
                  TVS2=int2str(TAVS2);
                  str2=TVS2;    
                  str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTM2x2',')',Semi];
                  cmd = sprintf(str)
                  eval(cmd)  
 #          disp('  Nodal forces at the node')
                  strSolmB=num2str(Solm);
                  strB1=str1;
                  NrB2=Nr1+1;
                  strB2=num2str(NrB2);
                  Arv123=Arv2;
          endfor 
 #           disp('   Nodal forces at the node')
                 strB=['B(',strB1, koolonV, strB2, Kom,Arv1,')','=','sSolmF(', Arv1, koolonV, Arv2, Kom, strSolmB,')',Semi];
                  cmd = sprintf(strB)
                 eval(cmd) 
                 Nr1=Nr1+2;  
                endif
       endif  ### (ElArvSolmes -2 ) == 0
 endif  ### if (wfiLaiend(Node,1) != 0)
  
endfor  #% for i=1:NSARV
#
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
vorrandeid=size(spA);
spRidu=vorrandeid(1,1);
disp('----- ')
vorrandeidE=size(spA);
disp(' ')
spA_rows=vorrandeidE(1,1)
spA_rowsE=spA_rows;
spA_cols=vorrandeidE(1,2)
nnzE=nnz(spA);
spA_nnz=sprintf('non_zero_elements_in_spA: %d', nnz(spA))
disp('  ')
disp('-----The equilibrium equations are inserted into spA  ----  ')
disp('  ')
equilibrium_equations_rows=spA_rowsE-spA_rowsD
non_zero_elements_in_equilibrium_equations= nnzE-nnzD
%nnzE=nnz(spA)
disp('  ')
vorrandeid=size(spA);
spRidu=vorrandeid(1,1); 
Nr1=spRidu+1;
ABB=ABB;
ABBdim=size(ABB) ;
#
#disp(' rida 662 ')
#spA
#spA=sparse(NNK,NNK)
#pause
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
disp(' Side conditions (hinges) ')
From_rows=sprintf('Side conditions begin from row: %d', Nr1)
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
#
#ABB
for i=1:ABBdim(1,1)
# 
    if ABB(i,8) == 1
      str1=num2str(Nr1); 
      Arv1=num2str(1);
      TAVS2=eval(sprintf('%3i',ABB(i,5)));  
      TVS2=int2str(TAVS2);
      str2=TVS2;    
      str=['spA=spSisestaArv(spA',Kom, str1,Kom, str2, Kom, Arv1,')',Semi];
# 
      cmd = sprintf(str)
      eval(cmd); 
      Nr1=Nr1+1  ;  
    endif
endfor
Nr1=Nr1-1;
#
disp('----- ')
vorrandeidN=size(spA);
#
disp(' ')
spA_rows=vorrandeidN(1,1)
spA_rowsN=spA_rows;
spA_cols=vorrandeidN(1,2)
#
nnzN=nnz(spA);
spA_nnz=sprintf('non_zero_elements_in_spA: %d', nnz(spA))
disp('  ')
disp('-----Side conditions are inserted into spA  ----  ')
disp('  ')
side_condition_rows=spA_rowsN-spA_rowsE
non_zero_elements_in_side_condition= nnzN-nnzE
#nnzN=nnz(spA)
disp('  ')
#%
#disp('----- ')
#
vorrandeid=size(spA);
spRidu=vorrandeid(1,1); 
Nr1=spRidu+1;   
#
  ABB;
  tsolm=tsolm;
  TSARV=size(tsolm);
  NTSARV=TSARV(1,1);
toereaktsioonid_algavad=Nr1;  
# 
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
disp(' Restrtictions on support displacements ')
From_rows=sprintf('Restrtictions on support displacements begin from row: %d', Nr1)
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
#disp(' rida 723 ')
#spA
#spA=sparse(NNK,NNK)
#pause
#
#ABB
#tsolm

TSARV=size(tsolm);
NTSARV=TSARV(1,1);
ABBdim=size(ABB);
ABBMoot=ABBdim(1,1);
ToeSidemeteArv=0;
ToeSiireRealNR=toereaktsioonid_algavad-1;

for i=1:NTSARV
tsolNr =tsolm(i,1);
  for j=2:3
    if (tsolm(i,j) == 1)
       ToeSidemeteArv=ToeSidemeteArv+1;
   endif    
  endfor
endfor
ToeSidemeteArv;
disp('   ')
ToeSidemed=zeros(ToeSidemeteArv,6); 
ToeSolmes=zeros(NTSARV,6); 

ToeSolmes(:,1:3)=tsolm(:,1:3); # 

for i=1:NTSARV
tsolNr =tsolm(i,1);
LeidsinVarda=0;
  for j=1:ABBMoot
  if (LeidsinVarda == 0)
    if (tsolNr-ABB(j,6)) == 0  
      ToeSolmes(i,4:6)=ABB(j,1:3); 
      LeidsinVarda=1;
     endif   
   endif     
  endfor
endfor
#%%%%%%%%%
#ToeSolmes
#%%%%%%%%
nmr=0;
  for j=1:NTSARV
    for k=2:3  
        KusOnYks=tsolm(j,k);
        if ((KusOnYks-1) == 0) 
##        if (KusOnYks)   
           nmr=nmr+1;
        ToeSiireRealNR = ToeSiireRealNR +1; 
        ToeSidemed(nmr,1) = ToeSiireRealNR;
        ToeSidemed(nmr,2) = ToeSolmes(j,4); 
        ToeSidemed(nmr,3) = k-1;
        ToeSidemed(nmr,4:5) = ToeSolmes(j,5:6);  
        ToeSidemed(nmr,6) = ToeSolmes(j,1);
       endif       
     endfor
  endfor
#%%%%%%%%%%%%%
ToeSidemed;
#disp(' rida 786 ')
#%%%%%%%%%%%%%
 for i=1:ToeSidemeteArv   
     Nr11=ToeSidemed(i,1);
      str1=num2str(Nr11); 
      Arv1=num2str(1);
 UWFi=ToeSidemed(i,3);
#
 switch (UWFi);
case{1}    
VarrasS=ToeSidemed(i,2);
#SpTUv=SpToeSiirdeUvektor(NSARV,NEARV,VarrasS,krdn,selem);
#SpTFiV=SpToeSiirdeFiVektorT(VarrasN)
SpTwV=SpToeSiirdeWvektorT(1);
TAVS=ToeSidemed(i,4);
               TVS2=int2str(TAVS);
               str2=TVS2;    
               
               str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTwV',')',Semi];
               cmd = sprintf(str)
               eval(cmd)
#
    strB1=str1;
    SolmB=ToeSidemed(i,6);  # OK oli w
    strSolmB=num2str(SolmB);
    ToeSidemed(1,1);
    Arv1x=num2str(1); 
    shift1=tSiire(1,1,SolmB);
    
    strShift1=num2str(shift1);
%    disp('  ')
     strTxt1=['Support shift  '  'tSiire(' Arv1x Kom Arv1 Kom strSolmB ')' ' = ' strShift1 '  in x direction at the node ' strSolmB ]; 
    disp(strTxt1)
      strB=['B(',strB1, Kom,Arv1,')','=','tSiire(', Arv1x, Kom, Arv1, Kom, strSolmB,')',Semi];
      cmd = sprintf(strB)
      eval(cmd);  
     disp('  ')     
case{2}
VarrasS=ToeSidemed(i,2);
#SpTWv=SpToeSiirdeWvektor(NSARV,NEARV,VarrasS,krdn,selem);
SpTFiV=SpToeSiirdeFiVektorT(1);
%TAVS=ToeSidemed(i,5);
TAVS=ToeSidemed(i,4);
               TVS2=int2str(TAVS);
               str2=TVS2;    
               
               str=['spA=spInsertBtoA(spA',Kom, str1,Kom, str2,Kom, 'SpTFiV',')',Semi];
               cmd = sprintf(str)
               eval(cmd);
#
    strB1=str1;
    SolmB=ToeSidemed(i,6);
    strSolmB=num2str(SolmB);
    ToeSidemed(1,1);
    Arv1z=num2str(2);
    shift2=tSiire(2,1,SolmB);
    strShift2=num2str(shift2);
%    disp('  ')
    strTxt2=['Support shift  '  'tSiire(' Arv1z Kom Arv1 Kom strSolmB ')' ' = ' strShift2 '  in z direction at the node ' strSolmB ]; 
##    disp(strTxt2)
      strB=['B(',strB1, Kom,Arv1,')','=','tSiire(', Arv1z, Kom, Arv1, Kom, strSolmB,')',Semi];
      cmd = sprintf(strB)
      eval(cmd);  
    disp('  ')                     
#%%%%%
endswitch  
# 
 endfor # for ToeSidemeteArv

#%%%%%%%%%%%%%%%
#
vorrandeid=size(spA);
spRidu=vorrandeid(1,1);   
#
disp('----- ')
vorrandeidR=size(spA);
#
disp(' ')
spA_rows=vorrandeidR(1,1)
spA_rowsR=spA_rows;
spA_cols=vorrandeidR(1,2)
#
nnzR=nnz(spA);
spA_nnz=sprintf('non_zero_elements_in_spA: %d', nnz(spA))
disp('  ')
disp('-----Restrtictions equations are inserted into spA  ----  ')
disp('  ')
restrtictions_equations_rows=spA_rowsR-spA_rowsN
non_zero_elements_in_restrtictions_equations= nnzR-nnzN
#
disp('  ')
#
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
#
spA_rank = sprank(spA)
#if spA_rank != NNK
if ~(spA_rank==NNK)
error(' spA is rank deficient. Rank must be NNK (NNK is the number of unknowns)!') 
end
#
spA;
#
B;
BsuurusS=size(B);
Bsuurus=BsuurusS(1,1);
#
#disp(' Nr    B ')
disp('  ')
#for i=1:Bsuurus
##  %9.3e %9.3f
#vorrandiNr=i;
#disp(sprintf('%3i %7.3f',vorrandiNr,B(i,1)))
#endfor
disp('  ')
# 
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
disp(' Solving system of sparse equations. ')
disp(' X=spA\B; ')
#Support_reactions=sprintf('are the last rows: %d', restrtictions_equations_rows)
Support_reactions=sprintf('Support reactions begin from X row: %d', toereaktsioonid_algavad)
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
#
#%%%%%%%%%%%%%%%%%%%%
X=spA\B; % lahend  võrrandisõteemist spA*X=B
#%%%%%%%%%%%%%%%%%%%%%
#
#disp(' Nr     X ')
#disp('  ')
#for i=1:Bsuurus
##  %9.3e   %9.3f
#vorrandiNr=i;
#disp(sprintf('%3i   %9.3e',vorrandiNr,X(i,1)))
#endfor
##
#disp('  ')
#Support_reactions=sprintf('Support reactions begin from X row: %d', toereaktsioonid_algavad)
#
disp('  ')
#%disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
#ToeSidemeteArv  
##selemjl(:,1:8) # OK
#%%selemjl
for i=1:NEARV
elemendiN(i,1)=i;
siireVardaA(i,1:2)=selemjl(i,5:6);
siireVardaL(i,1:2)=selemjl(i,1:2);
joudVardaA(i,1:2)=selemjl(i,7:8);
joudVardaL(i,1:2)=selemjl(i,3:4);
endfor
#
siireVardaA;
siireVardaL;
joudVardaA;
joudVardaL;
#%%%%%%%%%%%%
for i=1:NEARV
JrN(i,1)=i;
endfor
#
JrNT=JrN';
#
disp('=======================================================')
disp('     Support displacements/(Cz, CMy) at nodes  ') 
disp('   X_No,  Element_No,  w/fi_No,  w_No,  fi_No,  node_No ')
disp('-------------------------------------------------------')
for i=1:ToeSidemeteArv
disp(sprintf('     %2i     %2i         %2i        %3i    %3i     %3i  '  ,ToeSidemed(i,1),    ToeSidemed(i,2), ToeSidemed(i,3), ToeSidemed(i,4), ToeSidemed(i,5), ToeSidemed(i,6)))
endfor
disp('-------------------------------------------------------')
#
#The_displacement_DaF-s_numbers_at_beginning_of_the_element=sprintf('elementNr, displacement: %d',JrN(i,1) siireVardaA)
disp('==================================================================')
disp(' Displacements and force numbers at beginning of the element  ') 
disp('   Nr,    w,   fi     Q,   M  ')
disp('------------------------------------------------------------------')
for i=1:NEARV
disp(sprintf('  %2i     %2i   %2i     %2i   %2i  '  ,JrN(i,1),    siireVardaA(i,1), siireVardaA(i,2), joudVardaA(i,1), joudVardaA(i,2)))
endfor
%
disp('-------------------------------------------------------')
disp(' ')
disp('==================================================================')
disp(' Displacements and forces numbers at end of the element  ') 
disp('   Nr,    w,   fi      Q,   M  ')
disp('------------------------------------------------------------------')
for i=1:NEARV
disp(sprintf('  %2i     %2i   %2i      %2i   %2i    '  ,JrN(i,1),    siireVardaL(i,1), siireVardaL(i,2), joudVardaL(i,1), joudVardaL(i,2)))
endfor
#
disp('-------------------------------------------------------')
disp(' ')
#  
#
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
disp(' Initial parmeter vectors for elements displacements and forces  ')
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
disp(' ')
disp('---- Descaling multiplier for displacements = 1/baasi0 ------  ')
disp(' ')
#
for i=1:NEARV
#%siireVardaA(i,1:2)
AlgPar(i,1:4)=[X(siireVardaA(i,1),1)/baasi0  X(siireVardaA(i,2),1)/baasi0  ...
 X(joudVardaA(i,1),1)  X(joudVardaA(i,2),1) ];
endfor
#
disp('===============================================================')
disp(' Unscaled initial parmeter vector ')
disp('Element Nr   w          fi             Q          M ')
disp('---------------------------------------------------------------')
for i=1:NEARV
disp(sprintf('  %2i   %+9.3e   %+9.3e    %9.3f  %9.3f',JrN(i,1),  AlgPar(i,1),   AlgPar(i,2),   AlgPar(i,3),  AlgPar(i,4)))
#AlgPar
endfor
disp('---------------------------------------------------------------')
#
JJ=0;
for i=1:NEARV
#siireVardaL(i,1:2)
JJ=JJ+1;
SiireJoud(JJ,1:4)=[X(siireVardaA(i,1),1)/baasi0  X(siireVardaA(i,2),1)/baasi0  X(joudVardaA(i,1),1)  X(joudVardaA(i,2),1)];
JJ=JJ+1;
SiireJoud(JJ,1:4)=[X(siireVardaL(i,1),1)/baasi0  X(siireVardaL(i,2),1)/baasi0   X(joudVardaL(i,1),1)  X(joudVardaL(i,2),1) ];
endfor
#
NEARV2=NEARV*2;
Nn=0;
Nn2=0;
for i=1:NEARV
Nn=Nn+1;
Nn2=Nn2+1;
JrN2(Nn,1)=Nn2;
Nn=Nn+1;
JrN2(Nn,1)=Nn2;
endfor
#
#SiireJoud
#
#disp('=============================================================')
#disp('      Unscaled daf-s at beginning/end of elements ')
#disp('Element Nr   w          fi             Q          M ')
#disp('-------------------------------------------------------------')
#for i=1:NEARV2
#disp(sprintf('  %2i   %9.3e   %9.3e    %9.3f  %9.3f',JrN2(i,1),  SiireJoud(i,1),   SiireJoud(i,2),  #SiireJoud(i,3),  SiireJoud(i,4)))
#endfor
#disp('-------------------------------------------------------------')
disp('  ')
#disp('============================================================================')
#
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
##disp('Support reactions, toereaktsioonid_algavad  ')
Support_reactions=sprintf('Support reactions begin from row X: %d', toereaktsioonid_algavad)
disp('%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ')
disp('  ')
#ToeSidemed 6 ja 3
disp('==============================================================')
disp(' Nr         X          Node   (Cz <=> 1 / Cy(moment) <=> 2) ')
disp('--------------------------------------------------------------')
disp('  ')
for i=toereaktsioonid_algavad:Bsuurus
#  %9.3e   %9.3f
vorrandiNr=i;  %%9.3f
disp(sprintf('%3i   %+12.6e    %3i            %3i ',vorrandiNr,X(i,1),ToeSidemed(i-toereaktsioonid_algavad+1,6),ToeSidemed(i-toereaktsioonid_algavad+1,3) ))
endfor
#
disp('--------------------------------------------------------------')
disp('  ')

AlgParm=AlgPar;
#
endfunction
