Changes in explore/asymint.py [49eb251:1820208] in sasmodels
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
explore/asymint.py
r49eb251 r1820208 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( a*qa/2)89 siB = env.sas_sinx_x( b*qb/2)90 siC = env.sas_sinx_x( c*qc/2)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) 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, Fq96 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*c102 dra, drb, drc = slda-SLD_SOLVENT, sldb-SLD_SOLVENT, sldc-SLD_SOLVENT103 Aa, Ab, Ac = b*c, a*c, a*b104 Ta, Tb, Tc = a + 2*da, b + 2*db, c + 2*dc105 drVa, drVb, drVc = dra*a*Aa, drb*b*Ab, drc*c*Ac106 drVTa, drVTb, drVTc = dra*Ta*Aa, drb*Tb*Ab, drc*Tc*Ac107 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*siC115 + (drVTa*siAt-drVa*siA)*siB*siC116 + siA*(drVTb*siBt-drVb*siB)*siC117 + siA*siB*(drVTc*siCt-drVc*siC))118 Fq.__doc__ = "core-shell parallelepiped a=%g, b=%g c=%g"%(a, b, c)119 volume = a*b*c + 2*da*Aa + 2*db*Ab + 2*dc*Ac120 norm = 1/(volume*10000)121 95 return norm, Fq 122 96 … … 210 184 NORM, KERNEL = make_parallelepiped(A, B, C) 211 185 NORM_MP, KERNEL_MP = make_parallelepiped(A, B, C, env=MPenv) 212 elif shape == 'core_shell_parallelepiped':213 #A, B, C = 4450, 14000, 47214 #A, B, C = 445, 140, 47 # integer for the sake of mpf215 A, B, C = 6800, 114, 1380216 DA, DB, DC = 2300, 21, 58217 SLDA, SLDB, SLDC = "5", "-0.3", "11.5"218 if 1: # C shortest219 B, C = C, B220 DB, DC = DC, DB221 SLDB, SLDC = SLDC, SLDB222 elif 0: # C longest223 A, C = C, A224 DA, DC = DC, DA225 SLDA, SLDC = SLDC, SLDA226 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)228 186 elif shape == 'paracrystal': 229 187 LATTICE = 'bcc' … … 384 342 print("gauss-150", *gauss_quad_2d(Q, n=150)) 385 343 print("gauss-500", *gauss_quad_2d(Q, n=500)) 386 print("gauss-1025", *gauss_quad_2d(Q, n=1025))387 print("gauss-2049", *gauss_quad_2d(Q, n=2049))388 344 #gridded_2d(Q, n=2**8+1) 389 345 gridded_2d(Q, n=2**10+1) 390 #gridded_2d(Q, n=2**1 2+1)346 #gridded_2d(Q, n=2**13+1) 391 347 #gridded_2d(Q, n=2**15+1) 392 if shape not in ('paracrystal', 'core_shell_parallelepiped'): 393 # adaptive forms on models for which the calculations are fast enough 348 if shape != 'paracrystal': # adaptive forms are too slow! 394 349 print("dblquad", *scipy_dblquad_2d(Q)) 395 350 print("semi-romberg-100", *semi_romberg_2d(Q, n=100))
Note: See TracChangeset
for help on using the changeset viewer.