pulisco la memoria

kill(all);

definisco un dimensionamento di riferimento
per le valutazioni numeriche

dim : [
    dm  = 38,
    t   = 2 ,
    eta = 2 ,
    a   = 400 ,
    b   = 300 ,
    E   = 210000 ,
    nu  = 0.3
];
dim : append ( dim , ev([ 
    l = sqrt(a^2+b^2) ,
    A = %pi * dm * t ,
    J = integrate( (dm/2)^3 * cos(theta)^2 * t , theta, 0 , 2*%pi ),
    G = E / 2 / (1+nu)
] , dim) );

comando append: accoda due liste in un'unica complessiva

append( [ pippo, pluto ] , [ 1, 2, 3 ,paperino ] );

definisco delle funzioni ausiliarie
per descrivere l'andamento delle caratteristiche
di sollecitazione sui tratti

funzione N1 vale 1 sul primo estremo di un tratto
e 0 sul secondo estremo, con andamento lineare.
Definisco N1 come funzione di una ascissa
curvilinea adimensionale xi, 
xi=0 al primo estremo
xi=1 al secondo estremo

N1 : 1- xi;

funzione N2 vale 0 sul primo estremo di un tratto
e 1 sul secondo estremo, con andamento lineare.

N2 : xi;

se voglio definire un andamento lineare
di una funzione f da valore f1 sul primo estremo
a valore f2 sul secondo estremo

f1*N1+f2*N2;

inserisco i diagrammi di azione interna
della struttura

            |---P---|-----CB-----|-----CA-----|-----CO-----| 

N__OA: +P      +CB/b        +0           -CO/b        $

N__AB: -P*l/a  -CB*a/b/l    -CA*b/a/l    +CO*l/a/b    $

N__BO: +P*b/a  +0           +CA/a        -CO/a        $

Mf_OA: +0      +0           + 0*N1+CA*N2 +CO*N1+ 0*N2 $

Mf_AB: +0      + 0*N1+CB*N2 +CA*N1+ 0*N2 +0           $

Mf_BO: +0      +CB*N1+ 0*N2 +0           + 0*N1+CO*N2 $

T__OA: +0      +0           -CA/a        +CO/a        $

T__AB: +0      -CB/l        +CA/l        +0           $

T__BO: +0      +CB/b        +0           -CO/b        $

definiamo l'energia potenziale elastica sui tre tratti

noto che sul tratto OA,  dl = a * dxi

U_OA : integrate(
    (
            N__OA^2/2/E/A
        +   Mf_OA^2/2/E/J
        +   T__OA^2/2/G/A * eta
    ) * a
, xi , 0, 1
)$

U_AB : integrate(
    (
            N__AB^2/2/E/A
        +   Mf_AB^2/2/E/J
        +   T__AB^2/2/G/A * eta
    ) * l
, xi , 0, 1
)$

U_BO : integrate(
    (
            N__BO^2/2/E/A
        +   Mf_BO^2/2/E/J
        +   T__BO^2/2/G/A * eta
    ) * b
, xi , 0, 1
)$

calcolo l'energia potenziale totale

U : fullratsimp ( U_OA + U_AB + U_BO );

teorema di castigliano applicato alle reazioni vincolari 
interne

diff ( U , CB , 1 ) = thetaB

ove thetaB è la variazione di angolo al vertice B rispetto
all'indeformata

in coerenza con il vincolo di incastro interno, tale
variazione deve essere nulla

eqns : [
    diff(U,CB)=0,
    diff(U,CA)=0,
    diff(U,CO)=0
]$

risolvo il sistema di equazioni trovando CA,CB, CO

caso_iperstatico : linsolve (eqns, [CA,CB,CO] )$

definisco caso isostatico per confronto
le cerniere NON trasmettono coppie

caso_isostatico  : [ CA=0, CB=0 , CO=0];

ho due diverse definizioni dell'energia interna
per i due casi

U_iperstatico : ev(U, caso_iperstatico , infeval);

controllo che U_iperstatico NON sia dipendente da CA,
ossia sia una costante in CA,
ossia i termini CA *o* si autoeliminino *o* sono assenti

diff(U_iperstatico, CA);

diff(U_iperstatico, CB);

diff(U_iperstatico, CO);

calcolo l'energia interna nel caso isostatico

U_isostatico : ev(U, caso_isostatico , infeval);

utilizzando Castigliano, calcolo i cedimenti delta
nei due casi

delta_iper : diff(U_iperstatico, P);

delta_iso  : diff(U_isostatico,  P);