function [K,MC,Ge,Ce,V,DP,btotg,B,D,detJt]=Penta15(E,NU,rho,coordor,n27,Ns27,Nt27,Nu27,Ns21,Nt21,Nu21,GyroFlag)

D = (E/((1+NU)*(1-2*NU)))*[1-NU NU NU 0 0 0 ;
     NU 1-NU NU 0 0 0 ;
     NU NU 1- NU 0 0 0 ;
     0 0 0 (1-2*NU)/2 0 0 ;
     0 0 0 0 (1- 2*NU)/2 0 ;
     0 0 0 0 0 (1- 2*NU)/2];

% aggiungo coordinate nodi condensati
coord20=[coordor(1:3,:);coordor(3,:);coordor(4:6,:);coordor(6,:);coordor(7:8,:);...
        coordor(3,:);coordor(9:11,:);coordor(6,:);coordor(12:15,:);coordor(15,:);];
    
[MC20,V,DP,Ge20,Ce20]=Masses20(rho,coord20,n27,Ns27,Nt27,Nu27,GyroFlag);

[btotg, detJt,Tra, detra2, B]=Bmatrix_P3x7(coordor,Ns21,Nt21,Nu21);

ACell = repmat({Tra}, 1, 15);
Tr = blkdiag(ACell{:});

B=B/V;

Melc=reorder_dof20(MC20);

% condense the values of the degenerated nodes to the existing ones
[MC]=condense2ord(Melc);

if GyroFlag == 1
    [Ge]=condense2ord(Ge20);
    [Ce]=condense2ord(Ce20);
    Ge=Ge-Ge.';
else
    Ge=0;
    Ce=0;
end

[Ktot]= Kbrick_P3x7(btotg,D,detJt,detra2);
K=Tr'*Ktot*Tr;

end
