kill(all);

eq. tr (1)

define(
    srr_da_xy(t,sxx,syy,sxy) ,
    -trigreduce(+ cos(t)*(-syy*cos(t)-sxy*sin(t)) + sin(t)*(-cos(t)*sxy-sxx*sin(t)))
);

eq tt (1)

define(
    srt_da_xy(t,sxx,syy,sxy),
    -trigreduce(+ cos(t)*(+syy*sin(t)-sxy*cos(t)) + sin(t)*(+sin(t)*sxy-sxx*cos(t)))
);

eq tt (2)

define(
    stt_da_xy(t,sxx,syy,sxy),
    +trigreduce(+ cos(t)*(+sxx*cos(t)-sxy*sin(t)) + sin(t)*(-cos(t)*sxy+syy*sin(t)))
);

definizione far field = stato tensionale uniforme su sigma y = 1

stt_ff : stt_da_xy(t,0,1,0);
srr_ff : srr_da_xy(t,0,1,0);
srt_ff : srt_da_xy(t,0,1,0);

srr(phi) := 1/r * diff(phi,r,1)+1/r^2*diff(phi,t,2);
stt(phi) :=  diff(phi,r,2);
srt(phi) := -diff(1/r*diff(phi,t,1),r,1);

nabla2(phi):=diff(phi,r,2)+1/r*diff(phi,r,1)+1/r^2*diff(phi,t,2);

nabla4(phi):=(nabla2(nabla2(phi)));

termini phi da considerarsi

phiterms : append ( [r^2,log(r)] , ev( [r^(-n+2),r^(n),r^(-n)]*cos(n*t)  ,n=2));

unks   : makelist (A[i],         i,1,length(phiterms));
myphi  : sum(unks[i]*phiterms[i],i,1,length(phiterms));

fullratsimp(nabla4(myphi));

eq:[];

define(
    mysrr(r,t),
    srr(myphi)
);

define(
    mystt(r,t),
    stt(myphi)
);

define(
    mysrt(r,t),
    srt(myphi)
);

condizione di bordo interno scarico
e stato a remoto uniassiale trattivo

eqns:[
   mysrr(1,t)=0,
   mysrt(1,t)=0,
   limit(mysrr(r,t),r,inf)=srr_ff,
   limit(mysrt(r,t),r,inf)=srt_ff,
   limit(mystt(r,t),r,inf)=stt_ff
];

tali forme devono tuttavia valere per ogni t; 
i termini della matrice dei coefficienti (aumentata)
devono risultare tutti nulli

mat:augcoefmatrix(eqns,[cos(2*t),sin(2*t)]);

creo una equazione per ogni termine

neweqns:[]$
for i:1 thru 5 do(
    for j:1 thru 3 do(
        neweqns : append (neweqns,[mat[i,j]=0])
    )
)$
neweqns;

linsolve(neweqns,unks),globalsolve=true;

wxplot2d(mystt(1,t),[t,0,2*%pi]);