Changeset 49eb251 in sasmodels for explore/asymint.py
 Timestamp:
 Nov 6, 2017 10:49:46 AM (5 years ago)
 Branches:
 master, core_shell_microgels, magnetic_model, ticket1257vesicleproduct, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
 Children:
 74768cb, d8ac2ad
 Parents:
 393facf
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

explore/asymint.py
r1820208 r49eb251 86 86 a, b, c = env.mpf(a), env.mpf(b), env.mpf(c) 87 87 def Fq(qa, qb, qc): 88 siA = env.sas_sinx_x( 0.5*a*qa/2)89 siB = env.sas_sinx_x( 0.5*b*qb/2)90 siC = env.sas_sinx_x( 0.5*c*qc/2)88 siA = env.sas_sinx_x(a*qa/2) 89 siB = env.sas_sinx_x(b*qb/2) 90 siC = env.sas_sinx_x(c*qc/2) 91 91 return siA * siB * siC 92 92 Fq.__doc__ = "parallelepiped a=%g, b=%g c=%g"%(a, b, c) 93 93 volume = a*b*c 94 94 norm = CONTRAST**2*volume/10000 95 return norm, Fq 96 97 def make_core_shell_parallelepiped(a, b, c, da, db, dc, slda, sldb, sldc, env=NPenv): 98 a, b, c = env.mpf(a), env.mpf(b), env.mpf(c) 99 da, db, dc = env.mpf(da), env.mpf(db), env.mpf(dc) 100 slda, sldb, sldc = env.mpf(slda), env.mpf(sldb), env.mpf(sldc) 101 drV0 = CONTRAST*a*b*c 102 dra, drb, drc = sldaSLD_SOLVENT, sldbSLD_SOLVENT, sldcSLD_SOLVENT 103 Aa, Ab, Ac = b*c, a*c, a*b 104 Ta, Tb, Tc = a + 2*da, b + 2*db, c + 2*dc 105 drVa, drVb, drVc = dra*a*Aa, drb*b*Ab, drc*c*Ac 106 drVTa, drVTb, drVTc = dra*Ta*Aa, drb*Tb*Ab, drc*Tc*Ac 107 def Fq(qa, qb, qc): 108 siA = env.sas_sinx_x(a*qa/2) 109 siB = env.sas_sinx_x(b*qb/2) 110 siC = env.sas_sinx_x(c*qc/2) 111 siAt = env.sas_sinx_x(Ta*qa/2) 112 siBt = env.sas_sinx_x(Tb*qb/2) 113 siCt = env.sas_sinx_x(Tc*qc/2) 114 return (drV0*siA*siB*siC 115 + (drVTa*siAtdrVa*siA)*siB*siC 116 + siA*(drVTb*siBtdrVb*siB)*siC 117 + siA*siB*(drVTc*siCtdrVc*siC)) 118 Fq.__doc__ = "coreshell parallelepiped a=%g, b=%g c=%g"%(a, b, c) 119 volume = a*b*c + 2*da*Aa + 2*db*Ab + 2*dc*Ac 120 norm = 1/(volume*10000) 95 121 return norm, Fq 96 122 … … 184 210 NORM, KERNEL = make_parallelepiped(A, B, C) 185 211 NORM_MP, KERNEL_MP = make_parallelepiped(A, B, C, env=MPenv) 212 elif shape == 'core_shell_parallelepiped': 213 #A, B, C = 4450, 14000, 47 214 #A, B, C = 445, 140, 47 # integer for the sake of mpf 215 A, B, C = 6800, 114, 1380 216 DA, DB, DC = 2300, 21, 58 217 SLDA, SLDB, SLDC = "5", "0.3", "11.5" 218 if 1: # C shortest 219 B, C = C, B 220 DB, DC = DC, DB 221 SLDB, SLDC = SLDC, SLDB 222 elif 0: # C longest 223 A, C = C, A 224 DA, DC = DC, DA 225 SLDA, SLDC = SLDC, SLDA 226 NORM, KERNEL = make_core_shell_parallelepiped(A, B, C, DA, DB, DC, SLDA, SLDB, SLDC) 227 NORM_MP, KERNEL_MP = make_core_shell_parallelepiped(A, B, C, DA, DB, DC, SLDA, SLDB, SLDC, env=MPenv) 186 228 elif shape == 'paracrystal': 187 229 LATTICE = 'bcc' … … 342 384 print("gauss150", *gauss_quad_2d(Q, n=150)) 343 385 print("gauss500", *gauss_quad_2d(Q, n=500)) 386 print("gauss1025", *gauss_quad_2d(Q, n=1025)) 387 print("gauss2049", *gauss_quad_2d(Q, n=2049)) 344 388 #gridded_2d(Q, n=2**8+1) 345 389 gridded_2d(Q, n=2**10+1) 346 #gridded_2d(Q, n=2**1 3+1)390 #gridded_2d(Q, n=2**12+1) 347 391 #gridded_2d(Q, n=2**15+1) 348 if shape != 'paracrystal': # adaptive forms are too slow! 392 if shape not in ('paracrystal', 'core_shell_parallelepiped'): 393 # adaptive forms on models for which the calculations are fast enough 349 394 print("dblquad", *scipy_dblquad_2d(Q)) 350 395 print("semiromberg100", *semi_romberg_2d(Q, n=100))
Note: See TracChangeset
for help on using the changeset viewer.