-
Notifications
You must be signed in to change notification settings - Fork 0
/
nonlinear_moduli.py
47 lines (42 loc) · 23.7 KB
/
nonlinear_moduli.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
def P0(K0, mu0, K1, mu1, c=0):
return [[-(40824*c**5*(K0 - K1)**3*(K0 + 2*mu0)**3*(mu0 - mu1)**3 + 108*c**4*(K0 - K1)**2*(mu0 - mu1)**2*(4*K0*mu0**2*(2*mu0*(3049*mu0 - 1159*mu1) + 189*K1*(mu0 - 36*mu1)) + 8*mu0**3*(K1*(353*mu0 - 2243*mu1) + 1512*mu0*(mu0 - mu1)) + 567*K0**4*(3*mu0 + 2*mu1) - 189*K0**3*(3*K1*(mu0 + 4*mu1) - 4*mu0*(13*mu0 + 7*mu1)) + 756*K0**2*mu0*(mu0*(29*mu0 + 6*mu1) - 2*K1*(mu0 + 9*mu1))) + 54*c**3*(mu0 - mu1)*(567*K0**6*(3*mu0 + 2*mu1)**2 + 189*K0**5*(3*mu0 + 2*mu1)*(-9*K1*(mu0 + 4*mu1) + 2*mu0*(29*mu0 + 6*mu1)) + 32*mu0**3*(448*mu0**3*(mu0 - mu1)**2 + K1**3*(554*mu0**2 + 467*mu0*mu1 - 2596*mu1**2) + 2*K1**2*mu0*(403*mu0**2 + 1084*mu0*mu1 - 1487*mu1**2)) + 2*K0**2*mu0*(-240*K1*mu0**2*(247*mu0**2 + 451*mu0*mu1 - 383*mu1**2) + 189*K1**3*(57*mu0**2 - 64*mu0*mu1 - 168*mu1**2) + 160*mu0**3*(341*mu0**2 - 304*mu0*mu1 - 37*mu1**2) - 3024*K1**2*mu0*(2*mu0**2 - 19*mu0*mu1 - 33*mu1**2)) + 3*K0**3*(16*mu0**3*(2843*mu0**2 - 436*mu0*mu1 - 1357*mu1**2) - 126*K1**2*mu0*(47*mu0**2 - 244*mu0*mu1 - 328*mu1**2) + 63*K1**3*(19*mu0**2 - 38*mu0*mu1 - 56*mu1**2) - 504*K1*mu0**2*(57*mu0**2 + 131*mu0*mu1 + 12*mu1**2)) - 63*K0**4*(-16*mu0**2*(92*mu0**2 + 31*mu0*mu1 - 23*mu1**2) + 45*K1**2*(mu0**2 - 8*mu0*mu1 - 8*mu1**2) + 6*K1*mu0*(77*mu0**2 + 296*mu0*mu1 + 152*mu1**2)) + 8*K0*mu0**2*(6720*mu0**3*(mu0 - mu1)**2 - 16*K1*mu0**2*(403*mu0**2 + 1084*mu0*mu1 - 1487*mu1**2) + 63*K1**3*(73*mu0**2 - 21*mu0*mu1 - 252*mu1**2) + 18*K1**2*mu0*(227*mu0**2 + 596*mu0*mu1 + 227*mu1**2))) + c*(5103*K0**5*(3*K1 + 4*mu0)*(3*mu0 + 2*mu1)**3 - 10206*K0**4*(3*K1 + 4*mu0)*(3*mu0 + 2*mu1)**2*(-8*mu0**2 + 5*K1*mu1 - 2*mu0*mu1) - 567*K0**3*(3*K1 + 4*mu0)*(3*mu0 + 2*mu1)*(-12*K1*mu0*(9*mu0**2 - 73*mu0*mu1 - 86*mu1**2) - 16*mu0**2*(69*mu0**2 + 37*mu0*mu1 - 31*mu1**2) + 9*K1**2*(3*mu0**2 - 16*mu0*mu1 - 12*mu1**2)) + 64*mu0**3*(3600*K1*mu0**2*(mu0 - mu1)**3 - 27*K1**3*(286*mu0**3 + 402*mu0**2*mu1 - 1662*mu0*mu1**2 - 1651*mu1**3) - 54*K1**2*mu0*(317*mu0**3 + 729*mu0**2*mu1 - 1359*mu0*mu1**2 - 1437*mu1**3) + 32*mu0**3*(529*mu0**3 + 933*mu0**2*mu1 - 303*mu0*mu1**2 - 1159*mu1**3)) + 144*K0*mu0**2*(48*K1*mu0**2*(367*mu0**3 + 579*mu0**2*mu1 - 1209*mu0*mu1**2 - 1487*mu1**3) - 18*K1**2*mu0*(1089*mu0**3 + 513*mu0**2*mu1 - 2193*mu0*mu1**2 - 1159*mu1**3) + 160*mu0**3*(219*mu0**3 + 379*mu0**2*mu1 - 225*mu0*mu1**2 - 373*mu1**3) - 189*K1**3*(52*mu0**3 - 16*mu0**2*mu1 - 249*mu0*mu1**2 - 162*mu1**3)) - 756*K0**2*mu0*(3*K1 + 4*mu0)*(-12*K1*mu0*(9*mu0**3 - 127*mu0**2*mu1 - 398*mu0*mu1**2 - 234*mu1**3) + 9*K1**2*(48*mu0**3 - 79*mu0**2*mu1 - 236*mu0*mu1**2 - 108*mu1**3) + 16*mu0**2*(-191*mu0**3 - 247*mu0**2*mu1 + 142*mu0*mu1**2 + 171*mu1**3))) + (3*K1 + 4*mu0)**2*(1701*K0**4*(3*mu0 + 2*mu1)**3 + 378*K0**3*(3*mu0 + 2*mu1)**2*(-9*K1*mu0 + 42*mu0**2 - 6*K1*mu1 + 58*mu0*mu1) + 1512*K0**2*mu0*(-9*K1*mu0 + 18*mu0**2 - 6*K1*mu1 + 22*mu0*mu1)*(6*mu0**2 + 13*mu0*mu1 + 6*mu1**2) + 2016*K0*mu0**2*(2*mu0 + 3*mu1)**2*(10*mu0*(mu0 + mu1) - 3*K1*(3*mu0 + 2*mu1)) + 32*mu0**3*(4*mu0*(87*mu0**3 + 579*mu0**2*mu1 + 681*mu0*mu1**2 + 403*mu1**3) - 3*K1*(249*mu0**3 + 933*mu0**2*mu1 + 1587*mu0*mu1**2 + 731*mu1**3))) + 9*c**2*(567*K0**5*(3*K1*(3*mu0 - 8*mu1) + 4*mu0*(8*mu0 - 3*mu1))*(3*mu0 + 2*mu1)**2 + 567*K0**6*(3*mu0 + 2*mu1)**3 - 189*K0**4*(3*mu0 + 2*mu1)*(45*K1**2*(3*mu0**2 - 8*mu1**2) - 36*K1*mu0*(7*mu0**2 - 24*mu0*mu1 - 8*mu1**2) + 16*mu0**2*(-67*mu0**2 + 9*mu0*mu1 + 33*mu1**2)) + 21*K0**3*(27*K1**3*(63*mu0**3 + 96*mu0**2*mu1 - 156*mu0*mu1**2 - 128*mu1**3) + 64*mu0**3*(836*mu0**3 + 42*mu0**2*mu1 - 642*mu0*mu1**2 - 111*mu1**3) - 324*K1**2*mu0*(34*mu0**3 + 53*mu0**2*mu1 - 108*mu0*mu1**2 - 104*mu1**3) + 432*K1*mu0**2*(18*mu0**3 - 164*mu0**2*mu1 - 51*mu0*mu1**2 + 72*mu1**3)) + 16*K0*mu0**2*(160*mu0**3*(mu0 - mu1)**2*(299*mu0 + 373*mu1) + 189*K1**3*(18*mu0**3 + 226*mu0**2*mu1 - 81*mu0*mu1**2 - 288*mu1**3) - 18*K1**2*mu0*(359*mu0**3 - 2337*mu0**2*mu1 - 603*mu0*mu1**2 + 831*mu1**3) + 96*K1*mu0**2*(151*mu0**3 - 768*mu0**2*mu1 - 492*mu0*mu1**2 + 1109*mu1**3)) + 12*K0**2*mu0*(189*K1**3*(47*mu0**3 + 164*mu0**2*mu1 - 144*mu0*mu1**2 - 192*mu1**3) - 756*K1**2*mu0*(18*mu0**3 + 31*mu0**2*mu1 - 66*mu0*mu1**2 - 108*mu1**3) + 120*K1*mu0**2*(191*mu0**3 - 1665*mu0**2*mu1 - 393*mu0*mu1**2 + 1517*mu1**3) + 32*mu0**3*(3237*mu0**3 - 1661*mu0**2*mu1 - 3239*mu0*mu1**2 + 1663*mu1**3)) - 64*mu0**3*(200*K1*mu0**2*(mu0 - mu1)**3 - 1344*mu0**3*(mu0 - mu1)**2*(2*mu0 + 3*mu1) + 18*K1**2*mu0*(109*mu0**3 - 537*mu0**2*mu1 - 303*mu0*mu1**2 + 731*mu1**3) + 3*K1**3*(131*mu0**3 - 1863*mu0**2*mu1 - 342*mu0*mu1**2 + 2949*mu1**3))))/(7*(3*c*(K0 - K1) + 3*K1 + 4*mu0)**3*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3),
(-2*(-1 + c)*(8*mu0**3*(720*(-1 + c)*K1*mu0**2*(mu0 - mu1)**2*(69*mu0 + 91*mu1) - 320*mu0**3*(mu0 - mu1)**2*(69*mu0 + 91*mu1) + 27*K1**3*((103 + 3051*c + 1989*c**2 + 1857*c**3)*mu0**3 + (2251 + 3327*c - 807*c**2 - 4771*c**3)*mu0**2*mu1 + (-1 + c)**2*(3589 + 3971*c)*mu0*mu1**2 - 1057*(-1 + c)**3*mu1**3) + 36*K1**2*mu0*((-587 + 4086*c + 1989*c**2 + 1512*c**3)*mu0**3 - 3*(-907 - 874*c + 269*c**2 + 1512*c**3)*mu0**2*mu1 + 3*(-1 + c)**2*(1573 + 1512*c)*mu0*mu1**2 - 21*(-1 + c)**2*(-7 + 72*c)*mu1**3)) + 18*K0**2*mu0*(27*K1**3*((701 + 2307*c + 1683*c**2 + 559*c**3)*mu0**3 + (2207 + 1359*c - 1989*c**2 - 1577*c**3)*mu0**2*mu1 + 7*(-1 + c)**2*(269 + 211*c)*mu0*mu1**2 - 459*(-1 + c)**3*mu1**3) - 108*K1**2*mu0*((139 + 573*c + 537*c**2 + 501*c**3)*mu0**3 - 3*(-151 - 347*c + 97*c**2 + 401*c**3)*mu0**2*mu1 + 21*(-1 + c)**2*(37 + 43*c)*mu0*mu1**2 - 3*(-1 + c)**2*(-127 + 67*c)*mu1**3) - 12*K1*mu0**2*((4244 + 16068*c + 11187*c**2 + 3501*c**3)*mu0**3 - 3*(-3996 - 3852*c + 3547*c**2 + 4301*c**3)*mu0**2*mu1 + 3*(4564 - 5572*c - 4093*c**2 + 5101*c**3)*mu0*mu1**2 + 3*(1692 - 3636*c + 3911*c**2 - 1967*c**3)*mu1**3) + 16*mu0**3*(3*(76 + 12*c + 663*c**2 + 504*c**3)*mu0**3 + (2276 + 6372*c - 807*c**2 - 4536*c**3)*mu0**2*mu1 + (3164 - 252*c - 4353*c**2 + 4536*c**3)*mu0*mu1**2 + 3*(444 - 2052*c + 1057*c**2 - 504*c**3)*mu1**3)) + 72*K0*mu0**2*(-80*mu0**3*(mu0 - mu1)**2*((44 + 69*c)*mu0 + (36 + 91*c)*mu1) + 27*K1**3*((113 + 781*c + 549*c**2 + 307*c**3)*mu0**3 + (681 + 617*c - 477*c**2 - 821*c**3)*mu0**2*mu1 + 7*(-1 + c)**2*(107 + 103*c)*mu0*mu1**2 - 207*(-1 + c)**3*mu1**3) - 15*K1**2*mu0*((172 - 1257*c - 522*c**2 + 207*c**3)*mu0**3 - 3*(380 - 153*c - 274*c**2 + 47*c**3)*mu0**2*mu1 - 3*(-1 + c)**2*(252 + 113*c)*mu0*mu1**2 + 3*(-1 + c)**2*(108 + 91*c)*mu1**3) - 8*K1*mu0**2*((1438 + 2061*c + 1989*c**2 + 1512*c**3)*mu0**3 + 3*(282 + 1499*c - 269*c**2 - 1512*c**3)*mu0**2*mu1 + 3*(798 - 859*c - 1451*c**2 + 1512*c**3)*mu0*mu1**2 + 3*(774 - 1327*c + 1057*c**2 - 504*c**3)*mu1**3)) + 243*K0**5*(21*K1*((3 + 2*c)*mu0 - 2*(-1 + c)*mu1)**3 + 2*mu0*((378 + 756*c + 504*c**2 + 57*c**3)*mu0**3 + (756 + 252*c - 672*c**2 - 271*c**3)*mu0**2*mu1 + 7*(72 - 96*c - 24*c**2 + 53*c**3)*mu0*mu1**2 + (112 - 336*c + 336*c**2 - 157*c**3)*mu1**3)) + 27*K0**3*(189*K1**3*((3 + 2*c)*mu0 - 2*(-1 + c)*mu1)**3 - 18*K1**2*mu0*((1134 + 3693*c + 2862*c**2 + 1061*c**3)*mu0**3 - 3*(-1176 - 887*c + 1102*c**2 + 961*c**3)*mu0**2*mu1 + 21*(-1 + c)**2*(152 + 123*c)*mu0*mu1**2 - (-1 + c)**2*(-896 + 761*c)*mu1**3) + 8*mu0**3*(3*(672 + 2244*c + 2196*c**2 + 893*c**3)*mu0**3 + (7392 + 10524*c - 5724*c**2 - 8837*c**3)*mu0**2*mu1 + (8568 - 7644*c - 8316*c**2 + 9637*c**3)*mu0*mu1**2 + (3024 - 9612*c + 7452*c**2 - 3479*c**3)*mu1**3) + 72*K1*mu0**2*((-336 - 922*c - 573*c**2 + 81*c**3)*mu0**3 + (-812 - 94*c + 849*c**2 + 57*c**3)*mu0**2*mu1 - 7*(74 - 122*c - 3*c**2 + 51*c**3)*mu0*mu1**2 + 3*(-28 + 54*c - 99*c**2 + 73*c**3)*mu1**3)) - 486*K0**4*(21*K1**2*((3 + 2*c)*mu0 - 2*(-1 + c)*mu1)**3 - 4*mu0**2*((252 + 714*c + 561*c**2 + 113*c**3)*mu0**3 + (714 + 518*c - 663*c**2 - 439*c**3)*mu0**2*mu1 + 7*(88 - 104*c - 51*c**2 + 77*c**3)*mu0*mu1**2 + 3*(56 - 168*c + 153*c**2 - 71*c**3)*mu1**3) - 5*K1*mu0*(c**3*(89*mu0 - 29*mu1)*(mu0 - mu1)**2 + 14*mu1*(3*mu0 + 2*mu1)**2 + 3*c**2*(mu0 - mu1)**2*(45*mu0 + 47*mu1) + 42*c*(3*mu0**3 + 5*mu0**2*mu1 - 4*mu0*mu1**2 - 4*mu1**3)))))/(7*(3*c*(K0 - K1) + 3*K1 + 4*mu0)**3*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3),
((-1 + c)*(3*K0 + 4*mu0)*(81*K0**5*(3*(63 + 126*c + 147*c**2 + 41*c**3)*mu0**3 - 6*(-63 - 21*c - 7*c**2 + 44*c**3)*mu0**2*mu1 + 3*(84 - 112*c + 56*c**2 + 53*c**3)*mu0*mu1**2 + 2*(28 - 84*c + 112*c**2 - 9*c**3)*mu1**3) + 4*mu0**2*(-64*mu0**3*(mu0 - mu1)*((617 + 504*c + 252*c**2)*mu0**2 + 2*(503 + 126*c - 252*c**2)*mu0*mu1 + 9*(53 - 84*c + 28*c**2)*mu1**2) + 45*K1**3*((-415 + 741*c - 97*c**2 + 471*c**3)*mu0**3 - 3*(145 - 351*c - 153*c**2 + 359*c**3)*mu0**2*mu1 + 3*(-1 + c)**2*(5 + 247*c)*mu0*mu1**2 - 135*(-1 + c)**3*mu1**3) + 36*(-1 + c)*K1**2*mu0*((2075 + 669*c + 756*c**2)*mu0**3 + 3*(725 + 31*c - 756*c**2)*mu0**2*mu1 + 3*(-25 - 731*c + 756*c**2)*mu0*mu1**2 - 27*(25 - 53*c + 28*c**2)*mu1**3) - 48*K1*mu0**2*((2075 + 669*c + 756*c**2)*mu0**3 + 3*(725 + 31*c - 756*c**2)*mu0**2*mu1 + 3*(-25 - 731*c + 756*c**2)*mu0*mu1**2 - 27*(25 - 53*c + 28*c**2)*mu1**3)) + 3*K0*mu0*(45*K1**3*((-815 + 1353*c + 439*c**2 + 423*c**3)*mu0**3 - 3*(305 - 943*c + 271*c**2 + 367*c**3)*mu0**2*mu1 + 3*(-1 + c)**2*(25 + 311*c)*mu0*mu1**2 - 255*(-1 + c)**3*mu1**3) - 36*K1**2*mu0*((4075 + 393*c - 1811*c**2 + 843*c**3)*mu0**3 - 3*(-1525 + 1233*c + 9*c**2 + 283*c**3)*mu0**2*mu1 - 3*(-1 + c)**2*(125 + 277*c)*mu0*mu1**2 + 3*(-1 + c)**2*(-425 + 279*c)*mu1**3) - 64*mu0**3*(3*(291 + 589*c + 84*c**2)*mu0**3 - (379 + 1941*c + 756*c**2)*mu0**2*mu1 + (-2561 - 1419*c + 756*c**2)*mu0*mu1**2 + (-1433 + 1593*c - 252*c**2)*mu1**3) - 48*K1*mu0**2*((4075 + 6319*c + 2094*c**2 + 1512*c**3)*mu0**3 - 3*(-1525 - 681*c + 694*c**2 + 1512*c**3)*mu0**2*mu1 + 3*(-125 - 681*c - 706*c**2 + 1512*c**3)*mu0*mu1**2 + 3*(-425 + 227*c + 702*c**2 - 504*c**3)*mu1**3)) + 9*K0**2*(45*K1**3*((-100 + 174*c + 92*c**2 + 9*c**3)*mu0**3 - 6*(20 - 67*c + 39*c**2 + 8*c**3)*mu0**2*mu1 + 3*(-1 + c)**2*(5 + 23*c)*mu0*mu1**2 - 30*(-1 + c)**3*mu1**3) - 9*K1**2*mu0*((2000 + 4107*c - 214*c**2 + 1107*c**3)*mu0**3 - 3*(-800 - 1073*c + 1046*c**2 + 827*c**3)*mu0**2*mu1 + 3*(-1 + c)**2*(-100 + 547*c)*mu0*mu1**2 - 3*(-1 + c)**2*(200 + 89*c)*mu1**3) + 16*mu0**3*((584 + 321*c + 585*c**2 + 756*c**3)*mu0**3 + (3428 + 4357*c + 345*c**2 - 2268*c**3)*mu0**2*mu1 + (4692 + 823*c - 2445*c**2 + 2268*c**3)*mu0*mu1**2 + (1796 - 2001*c + 1515*c**2 - 756*c**3)*mu1**3) - 12*K1*mu0**2*((2000 + 9174*c + 5657*c**2 + 669*c**3)*mu0**3 + 3*(800 + 2166*c - 1737*c**2 - 1229*c**3)*mu0**2*mu1 + 3*(-100 - 906*c - 783*c**2 + 1789*c**3)*mu0*mu1**2 - 3*(200 + 818*c - 1801*c**2 + 783*c**3)*mu1**3)) - 27*K0**4*(-28*mu0*(3*mu0 + 2*mu1)**2*(12*mu0 + 13*mu1) + 3*c**3*(mu0 - mu1)**2*(3*K1*(67*mu0 + 38*mu1) + mu0*(-303*mu0 + 583*mu1)) - 21*c*(3*mu0 + 2*mu1)*(3*K1*(3*mu0 + 2*mu1)**2 + 4*mu0*(38*mu0**2 + 9*mu0*mu1 - 22*mu1**2)) + c**2*(-7272*mu0**4 + 5436*mu0**3*mu1 + 444*mu0**2*mu1**2 - 5608*mu0*mu1**3 + 9*K1*(122*mu0**3 + 474*mu0**2*mu1 + 261*mu0*mu1**2 + 18*mu1**3))) + 9*K0**3*(112*mu0**2*(171*mu0**3 + 492*mu0**2*mu1 + 453*mu0*mu1**2 + 134*mu1**3) + 3*c**3*(mu0 - mu1)**2*(4*mu0**2*(727*mu0 - 1287*mu1) + 9*K1**2*(11*mu0 + 94*mu1) + 9*K1*mu0*(33*mu0 + 247*mu1)) - 3*c*(12*K1*mu0*(811*mu0**3 + 507*mu0**2*mu1 - 822*mu0*mu1**2 - 496*mu1**3) + 9*K1**2*(689*mu0**3 + 978*mu0**2*mu1 + 177*mu0*mu1**2 - 94*mu1**3) + 16*mu0**2*(-1026*mu0**3 - 1437*mu0**2*mu1 + 177*mu0*mu1**2 + 536*mu1**3)) - c**2*(27*K1**2*(197*mu0**3 - 906*mu0**2*mu1 - 354*mu0*mu1**2 + 188*mu1**3) + 45*K1*mu0*(771*mu0**3 - 297*mu0**2*mu1 + 381*mu0*mu1**2 + 545*mu1**3) - 4*mu0**2*(8271*mu0**3 - 6753*mu0**2*mu1 - 7107*mu0*mu1**2 + 9089*mu1**3)))))/(21*(3*c*(K0 - K1) + 3*K1 + 4*mu0)**3*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3)],
[(-120*(-1 + c)*mu0**2*(mu0 - mu1)**2*(6*K0*mu0*((42 + 23*c)*mu0 + (28 - 23*c)*mu1) + 63*K0*K1*((3 + 2*c)*mu0 - 2*(-1 + c)*mu1) + 2*mu0*(4*mu0*((23 + 42*c)*mu0 + (47 - 42*c)*mu1) + 3*K1*((23 + 47*c)*mu0 - 47*(-1 + c)*mu1))))/(7*(3*c*(K0 - K1) + 3*K1 + 4*mu0)*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3),
(-567*K0**4*((3 + 2*c)*mu0 - 2*(-1 + c)*mu1)**3 + 8*mu0**3*(3*(-1 + c)*K1*(mu0 - mu1)*((-223 + 3261*c + 1512*c**2)*mu0**2 + 2*(2278 + 2559*c - 1512*c**2)*mu0*mu1 + 63*(109 - 133*c + 24*c**2)*mu1**2) - 20*mu0*((17 + 879*c + 504*c**2)*mu0**3 + 3*(383 + 121*c - 504*c**2)*mu0**2*mu1 + 3*(337 - 841*c + 504*c**2)*mu0*mu1**2 - 21*(37 - 61*c + 24*c**2)*mu1**3)) + 12*K0*mu0**2*(3*(-1 + c)*K1*(mu0 - mu1)*((907 + 4756*c + 1512*c**2)*mu0**2 + 3*(2447 + 836*c - 1008*c**2)*mu0*mu1 + 8*(719 - 908*c + 189*c**2)*mu1**2) - 2*mu0*((3368 + 13811*c + 9309*c**2 + 1512*c**3)*mu0**3 + (18986 + 1837*c - 16287*c**2 - 4536*c**3)*mu0**2*mu1 + 3*(3108 - 6169*c + 1549*c**2 + 1512*c**3)*mu0*mu1**2 + 3*(-1226 + 953*c + 777*c**2 - 504*c**3)*mu1**3)) + 54*K0**3*(21*(-1 + c)*K1*(mu0 - mu1)*((3 + 2*c)**2*mu0**2 + (27 + 6*c - 8*c**2)*mu0*mu1 + 2*(7 - 9*c + 2*c**2)*mu1**2) - mu0*((1176 + 2887*c + 2433*c**2 + 504*c**3)*mu0**3 + (3262 + 1609*c - 3359*c**2 - 1512*c**3)*mu0**2*mu1 + 7*(314 - 447*c - 83*c**2 + 216*c**3)*mu0*mu1**2 + (364 - 1367*c + 1507*c**2 - 504*c**3)*mu1**3)) + 18*K0**2*mu0*(3*(-1 + c)*K1*(mu0 - mu1)*(3*(323 + 559*c + 168*c**2)*mu0**2 + 2*(1436 + 293*c - 504*c**2)*mu0*mu1 + (1759 - 2263*c + 504*c**2)*mu1**2) - 2*mu0*((3653 + 8111*c + 7724*c**2 + 1512*c**3)*mu0**3 - 2*(-5193 - 2651*c + 5576*c**2 + 2268*c**3)*mu0**2*mu1 + 7*(967 - 1491*c - 124*c**2 + 648*c**3)*mu0*mu1**2 - 24*(-8 + 124*c - 179*c**2 + 63*c**3)*mu1**3)))/(7*(3*c*(K0 - K1) + 3*K1 + 4*mu0)*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3),
-((-1 + c)*(3*K0 + 4*mu0)*(27*K0**3*((420 + 221*c)*mu0**3 + (490 - 278*c)*mu0**2*mu1 + (35 - 107*c)*mu0*mu1**2 + 2*(-35 + 82*c)*mu1**3) - 4*mu0**2*(4*mu0*(mu0 - mu1)*((1235 + 84*c)*mu0**2 - 2*(-515 + 84*c)*mu0*mu1 + 3*(-55 + 28*c)*mu1**2) + 15*K1*(3*(157 + 39*c)*mu0**3 + (407 - 183*c)*mu0**2*mu1 + (-127 + 15*c)*mu0*mu1**2 + 51*(-1 + c)*mu1**3)) - 3*K0*mu0*(4*mu0*((855 + 143*c)*mu0**3 - 27*(95 + 47*c)*mu0**2*mu1 + (-2335 + 2109*c)*mu0*mu1**2 + (545 - 983*c)*mu1**3) + 15*K1*((955 + 389*c)*mu0**3 + (775 - 663*c)*mu0**2*mu1 + (-215 + 159*c)*mu0*mu1**2 + 115*(-1 + c)*mu1**3)) - 9*K0**2*(15*K1*((121 + 61*c)*mu0**3 + 2*(46 - 53*c)*mu0**2*mu1 + (-22 + 29*c)*mu0*mu1**2 + 16*(-1 + c)*mu1**3) - mu0*((2060 + 993*c)*mu0**3 + 27*(160 - 17*c)*mu0**2*mu1 + (1280 - 2061*c)*mu0*mu1**2 + 3*(-220 + 509*c)*mu1**3))))/(42*(3*c*(K0 - K1) + 3*K1 + 4*mu0)*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3)],
[(7200*(-1 + c)*mu0**3*(mu0 - mu1)**3)/(7*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3),
(-180*(-1 + c)*mu0*(mu0 - mu1)**2*(9*K0**2*(11*mu0 + 9*mu1) + 36*K0*mu0*(11*mu0 + 9*mu1) + 4*mu0**2*(69*mu0 + 91*mu1)))/(7*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3),
-(27*K0**2*mu0*((2025 + 1139*c + 336*c**2)*mu0**3 - 3*(-775 + 439*c + 336*c**2)*mu0**2*mu1 + 9*(-25 - 87*c + 112*c**2)*mu0*mu1**2 + (-625 + 961*c - 336*c**2)*mu1**3) + 72*K0*mu0**2*((1025 + 473*c + 252*c**2)*mu0**3 - 9*(-125 + 41*c + 84*c**2)*mu0**2*mu1 + 3*(-25 - 227*c + 252*c**2)*mu0*mu1**2 + (-325 + 577*c - 252*c**2)*mu1**3) + 27*K0**3*((500 + 319*c + 56*c**2)*mu0**3 - 24*(-25 + 18*c + 7*c**2)*mu0**2*mu1 + 3*(-25 - 31*c + 56*c**2)*mu0*mu1**2 - 2*(75 - 103*c + 28*c**2)*mu1**3) + 16*mu0**3*((2075 + 669*c + 756*c**2)*mu0**3 + 3*(725 + 31*c - 756*c**2)*mu0**2*mu1 + 3*(-25 - 731*c + 756*c**2)*mu0*mu1**2 - 27*(25 - 53*c + 28*c**2)*mu1**3))/(7*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3)]]
def P1(K0, mu0, K1, mu1, c=0):
return [[(3*K0 + 4*mu0)**3/(3*c*(K0 - K1) + 3*K1 + 4*mu0)**3,
0,
-(((-1 + c)*(3*K0 + 4*mu0)**3*(K0*(3*mu0 + 2*mu1) + mu0*(-5*K1 - 4*mu0 + 4*mu1))*(3*K0**2*(3*(3 + 9*c + 13*c**2)*mu0**2 + 6*(2 + c - 3*c**2)*mu0*mu1 + 4*(-1 + c)**2*mu1**2) + mu0**2*(75*(-1 + c)**2*K1**2 - 60*(-1 + c)*K1*((4 + c)*mu0 + mu1 - c*mu1) + 16*((13 + 9*c + 3*c**2)*mu0**2 - 3*(-3 + c + 2*c**2)*mu0*mu1 + 3*(-1 + c)**2*mu1**2)) + 3*K0*mu0*(-5*(-1 + c)*K1*(3*(1 + 4*c)*mu0 - 2*(-1 + c)*mu1) + 4*((9 + 32*c + 9*c**2)*mu0**2 + (12 + c - 13*c**2)*mu0*mu1 + 4*(-1 + c)**2*mu1**2))))/((3*c*(K0 - K1) + 3*K1 + 4*mu0)**3*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3))],
[0,
(25*mu0**2*(3*K0 + 4*mu0)**3)/((3*c*(K0 - K1) + 3*K1 + 4*mu0)*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**2),
(25*(-1 + c)*mu0**2*(3*K0 + 4*mu0)**3*(K0*(3*mu0 + 2*mu1) + mu0*(-5*K1 - 4*mu0 + 4*mu1)))/(2*(3*c*(K0 - K1) + 3*K1 + 4*mu0)*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3)],
[0,
0,
(125*mu0**3*(3*K0 + 4*mu0)**3)/(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3]]
def Q(K0, mu0, K1, mu1, c=0):
return [-((-1 + c)*(16*mu0**4*(108*(-1 + c)**2*K1**2*mu0*(mu0 - mu1)**2*((-335 + 1008*c)*mu0 - (2465 + 1008*c)*mu1) - 144*(-1 + c)*K1*mu0**2*(mu0 - mu1)**2*((-335 + 1008*c)*mu0 - (2465 + 1008*c)*mu1) + 64*mu0**3*(mu0 - mu1)**2*((-335 + 1008*c)*mu0 - (2465 + 1008*c)*mu1) + 27*K1**3*((1009 + 8957*c + 9075*c**2 + 13943*c**3 + 2016*c**4)*mu0**3 + (4253 + 18489*c + 7215*c**2 - 23909*c**3 - 6048*c**4)*mu0**2*mu1 + (-1 + c)**2*(13667 + 18085*c + 6048*c**2)*mu0*mu1**2 - (-1 + c)**3*(2071 + 2016*c)*mu1**3)) + 48*K0*mu0**3*(-72*(-1 + c)*K1*mu0**2*(mu0 - mu1)**2*((-140 + 1093*c + 1008*c**2)*mu0 - (2860 + 3893*c + 1008*c**2)*mu1) + 16*mu0**3*(mu0 - mu1)**2*((-280 + 1851*c + 3024*c**2)*mu0 - (5720 + 10251*c + 3024*c**2)*mu1) + 54*K1**3*(7*(151 + 750*c + 879*c**2 + 788*c**3 + 57*c**4)*mu0**3 - 3*(-1213 - 3130*c + 148*c**2 + 3796*c**3 + 399*c**4)*mu0**2*mu1 + 9*(-1 + c)**2*(764 + 958*c + 133*c**2)*mu0*mu1**2 - (-1 + c)**3*(1553 + 399*c)*mu1**3) - 27*K1**2*mu0*((728 - 169*c + 11370*c**2 + 7055*c**3 + 2016*c**4)*mu0**3 + (7176 + 8107*c - 5990*c**2 - 3245*c**3 - 6048*c**4)*mu0**2*mu1 + (-1 + c)**2*(-8136 - 2579*c + 6048*c**2)*mu0*mu1**2 - (-1 + c)**2*(-7232 - 6833*c + 2016*c**2)*mu1**3)) + 81*K0**4*(567*K1**3*((3 + 2*c)*mu0 - 2*(-1 + c)*mu1)**3 - 54*K1**2*mu0*(7*(243 + 720*c + 660*c**2 + 248*c**3 + 4*c**4)*mu0**3 - 2*(-2331 - 2095*c + 2300*c**2 + 2084*c**3 + 42*c**4)*mu0**2*mu1 + 4*(-1 + c)**2*(987 + 824*c + 21*c**2)*mu0*mu1**2 - 4*(-1 + c)**2*(-266 + 188*c + 7*c**2)*mu1**3) - 32*mu0**3*((-756 - 2268*c - 1632*c**2 + 119*c**3 + 399*c**4)*mu0**3 - 3*(924 + 402*c - 2687*c**2 - 1601*c**3 + 399*c**4)*mu0**2*mu1 + 3*(-1071 + 1322*c - 942*c**2 - 3321*c**3 + 399*c**4)*mu0*mu1**2 + (-1134 + 4758*c - 3603*c**2 + 5041*c**3 - 399*c**4)*mu1**3) + 72*K1*mu0**2*((126 + 1533*c + 3533*c**2 + 2487*c**3 + 196*c**4)*mu0**3 + (777 + 3401*c - 89*c**2 - 3501*c**3 - 588*c**4)*mu0**2*mu1 + (1218 + 1724*c - 3071*c**2 - 459*c**3 + 588*c**4)*mu0*mu1**2 + (504 - 1408*c - 373*c**2 + 1473*c**3 - 196*c**4)*mu1**3)) - 486*K0**5*(189*K1**2*((3 + 2*c)*mu0 - 2*(-1 + c)*mu1)**3 + 4*mu0**2*((-1134 - 2457*c - 630*c**2 + 835*c**3 + 196*c**4)*mu0**3 - 3*(1071 + 63*c - 1810*c**2 - 485*c**3 + 196*c**4)*mu0**2*mu1 + 3*(-924 + 1708*c - 540*c**2 - 1805*c**3 + 196*c**4)*mu0*mu1**2 + (-756 + 2772*c - 3180*c**2 + 3125*c**3 - 196*c**4)*mu1**3) - 9*K1*mu0*(7*(81 + 306*c + 336*c**2 + 148*c**3 + 4*c**4)*mu0**3 - 4*(-441 - 546*c + 449*c**2 + 517*c**3 + 21*c**4)*mu0**2*mu1 + 4*(399 - 336*c - 341*c**2 + 257*c**3 + 21*c**4)*mu0*mu1**2 + 4*(112 - 308*c + 202*c**2 + c**3 - 7*c**4)*mu1**3)) + 729*K0**6*(63*K1*((3 + 2*c)*mu0 - 2*(-1 + c)*mu1)**3 - 2*mu0*(7*(-81 - 108*c + 48*c**3 + 4*c**4)*mu0**3 + 2*(-567 + 189*c + 630*c**2 + 16*c**3 - 42*c**4)*mu0**2*mu1 + 4*(-189 + 378*c - 105*c**2 - 268*c**3 + 21*c**4)*mu0*mu1**2 - 4*(42 - 154*c + 210*c**2 - 176*c**3 + 7*c**4)*mu1**3)) + 54*K0**3*mu0*(27*K1**3*(7*(405 + 1146*c + 972*c**2 + 348*c**3 + 4*c**4)*mu0**3 - 4*(-1840 - 1536*c + 1788*c**2 + 1567*c**3 + 21*c**4)*mu0**2*mu1 + 4*(-1 + c)**2*(1675 + 1349*c + 21*c**2)*mu0*mu1**2 - 4*(-1 + c)**3*(370 + 7*c)*mu1**3) - 108*K1**2*mu0*((1386 + 5434*c + 7220*c**2 + 4139*c**3 + 196*c**4)*mu0**3 + (4967 + 8258*c - 4180*c**2 - 8457*c**3 - 588*c**4)*mu0**2*mu1 + (-1 + c)**2*(4808 + 5673*c + 588*c**2)*mu0*mu1**2 - (-1 + c)**2*(-1964 + 571*c + 196*c**2)*mu1**3) - 8*mu0**3*((-1344 - 5944*c - 11592*c**2 - 6721*c**3 + 2016*c**4)*mu0**3 + (-4448 + 8712*c + 50256*c**2 + 38083*c**3 - 6048*c**4)*mu0**2*mu1 + (-12272 - 20592*c - 40536*c**2 - 56003*c**3 + 6048*c**4)*mu0*mu1**2 + (-2936 + 31824*c + 1872*c**2 + 24641*c**3 - 2016*c**4)*mu1**3) + 144*K1*mu0**2*((-140 + 178*c + 2020*c**2 + 1918*c**3 + 399*c**4)*mu0**3 + (-620 + 356*c + 2055*c**2 - 594*c**3 - 1197*c**4)*mu0**2*mu1 + (85 + 5054*c - 1770*c**2 - 4566*c**3 + 1197*c**4)*mu0*mu1**2 - (200 + 338*c + 2305*c**2 - 3242*c**3 + 399*c**4)*mu1**3)) + 72*K0**2*mu0**2*(8*mu0**3*(mu0 - mu1)**2*((40 + 1624*c + 7563*c**2 + 3024*c**3)*mu0 - (4840 + 19624*c + 15963*c**2 + 3024*c**3)*mu1) + 27*K1**3*((2651 + 9806*c + 10431*c**2 + 5791*c**3 + 196*c**4)*mu0**3 - 3*(-2714 - 4234*c + 2281*c**2 + 4471*c**3 + 196*c**4)*mu0**2*mu1 + 3*(-1 + c)**2*(3471 + 3543*c + 196*c**2)*mu0*mu1**2 - (-1 + c)**3*(2419 + 196*c)*mu1**3) - 108*K1**2*mu0*((471 + 2173*c + 4615*c**2 + 3717*c**3 + 399*c**4)*mu0**3 - 3*(-849 - 1807*c + 260*c**2 + 1997*c**3 + 399*c**4)*mu0**2*mu1 + 3*(-1 + c)**2*(456 + 1075*c + 399*c**2)*mu0*mu1**2 - (-1 + c)**2*(-1739 - 645*c + 399*c**2)*mu1**3) + 18*K1*mu0**2*((-408 + 744*c + 4481*c**2 + 167*c**3 + 2016*c**4)*mu0**3 + (-6936 - 12792*c + 8357*c**2 + 17419*c**3 - 6048*c**4)*mu0**2*mu1 + (6696 + 27552*c - 4957*c**2 - 35339*c**3 + 6048*c**4)*mu0*mu1**2 - (6352 + 1504*c + 7881*c**2 - 17753*c**3 + 2016*c**4)*mu1**3))))/(42*(3*c*(K0 - K1) + 3*K1 + 4*mu0)**3*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3),
(-2*(-1 + c)*(mu0 - mu1)*(27*K0**4*(mu0 - mu1)*(105*K1*((3 + 2*c)*mu0 - 2*(-1 + c)*mu1) + mu0*(-7*(-21 + 16*c + 29*c**2)*mu0 + (98 - 88*c + 203*c**2)*mu1)) + 9*K0**3*mu0*(3*K1*(mu0 - mu1)*(7*(199 + 187*c + 29*c**2)*mu0 + (1312 - 1109*c - 203*c**2)*mu1) + mu0*((1946 - 1303*c - 3164*c**2)*mu0**2 + 2*(769 + 173*c + 3164*c**2)*mu0*mu1 + (-1384 + 957*c - 3164*c**2)*mu1**2)) + 2*mu0**4*(-4*mu0*(mu0 - mu1)*((-335 + 1008*c)*mu0 - (2465 + 1008*c)*mu1) + K1*((6409 + 10797*c + 13104*c**2)*mu0**2 - 2*(-1186 + 11637*c + 13104*c**2)*mu0*mu1 + 3*(-8527 + 4159*c + 4368*c**2)*mu1**2)) + 3*K0**2*mu0**2*(3*K1*((7069 + 8947*c + 3164*c**2)*mu0**2 - 2*(-341 + 7817*c + 3164*c**2)*mu0*mu1 + (-9851 + 6687*c + 3164*c**2)*mu1**2) + mu0*((9271 - 4242*c - 16044*c**2)*mu0**2 + (16063 + 3484*c + 32088*c**2)*mu0*mu1 - 2*(4267 - 379*c + 8022*c**2)*mu1**2)) + K0*mu0**3*(-2*mu0*((-7756 + 3909*c + 13104*c**2)*mu0**2 - 2*(13249 + 4749*c + 13104*c**2)*mu0*mu1 + 3*(5818 + 1863*c + 4368*c**2)*mu1**2) + 3*K1*(7*(2207 + 3586*c + 2292*c**2)*mu0**2 - 3*(-1299 + 15068*c + 10696*c**2)*mu0*mu1 + 2*(-18073 + 10051*c + 8022*c**2)*mu1**2))))/(7*(3*c*(K0 - K1) + 3*K1 + 4*mu0)*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3),
(12*(-1 + c)*mu0*(mu0 - mu1)**2*(63*c*(K0 + 2*mu0)**2*(9*K0 + 8*mu0)*(mu0 - mu1) - 5*(-9*K0**2*mu0*(mu0 - 121*mu1) + 270*K0**3*mu1 + 12*K0*mu0**2*(7*mu0 + 143*mu1) + 2*mu0**3*(67*mu0 + 493*mu1))))/(7*(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))**3)]
def relative_moduli(K0, mu0, l0, m0, n0, K1, mu1, l1, m1, n1, c=0):
K_rel = -(((K0 - K1)*(3*K0 + 4*mu0))/(3*c*(K0 - K1) + 3*K1 + 4*mu0))
mu_rel = (-5*mu0*(3*K0 + 4*mu0)*(mu0 - mu1))/(K0*((9 + 6*c)*mu0 - 6*(-1 + c)*mu1) + 4*mu0*((2 + 3*c)*mu0 - 3*(-1 + c)*mu1))
P0l, P0m, P0n = P0(K0, mu0, K1, mu1, c)
P1l, P1m, P1n = P1(K0, mu0, K1, mu1, c)
lg, mg, ng = Q(K0, mu0, K1, mu1)
l_rel = P0l[0]*l0 + P0l[1]*m0 + P0l[2]*n0 + P1l[0]*l1 + P1l[1]*m1 + P1l[2]*n1 + lg
m_rel = P0m[0]*l0 + P0m[1]*m0 + P0m[2]*n0 + P1m[0]*l1 + P1m[1]*m1 + P1m[2]*n1 + mg
n_rel = P0n[0]*l0 + P0n[1]*m0 + P0n[2]*n0 + P1n[0]*l1 + P1n[1]*m1 + P1n[2]*n1 + ng
return K_rel, mu_rel, l_rel, m_rel, n_rel
def effective_moduli(K0, mu0, l0, m0, n0, K1, mu1, l1, m1, n1, c):
K_rel, mu_rel, l_rel, m_rel, n_rel = relative_moduli(K0, mu0, l0, m0, n0, K1, mu1, l1, m1, n1, c)
K_eff = K0 + c*K_rel
mu_eff = mu0 + c*mu_rel
l_eff = l0 + c*l_rel
m_eff = m0 + c*m_rel
n_eff = n0 + c*n_rel
return K_eff, mu_eff, l_eff, m_eff, n_eff