ALPS 2 Tutorial:MC-09 Directed Worm Algorithm:Intensive Example 4
Contents
Waist correction
Starting the business
Bandstructure
Headers:
import numpy; import pyalps; import pyalps.dwa;
Try:
V0 = numpy.array([8.35, 8.35, 8.35]); wlen = numpy.array([843., 843., 843.]); a = 101 m = 86.99 L = 160 band = pyalps.dwa.bandstructure(V0, wlen, a, m, L);
See:
band
Consider a parabolic trapping of 10.5 Hz:
omega = numpy.array([10.5, 10.5, 10.5]); VT = pyalps.dwa.trappingFrequency(omega=omega, mass=m, wavelength=wlen) / band.t();
Parameters (no waist)
Parameters:
tof_phase = pyalps.dwa.tofPhase(time_of_flight=15.5, wavelength=[843.,843.,843.], mass=86.99) params=[] for seed in range(10,210,10): params.append( { 'LATTICE' : 'inhomogeneous simple cubic lattice' , 'MODEL' : 'boson Hubbard' , 'L' : 160 , 'Nmax' : 20 , 't' : 1. , 'U' : 8.1 , 'T' : 1. , 'mu_homogeneous' : 2. # Attempt to fix at N=280000 , lower bound: 2 , upper bound: 3 , 'mu' : 'mu_homogeneous - (0.00369*(x-(L-1)/2.)*(x-(L-1)/2.) + 0.00369*(y-(L-1)/2.)*(y-(L-1)/2.) + 0.00369*(z-(L-1)/2.)*(z-(L-1)/2.))' , 'tof_phase' : str(tof_phase) #, 'V0' : (8.35 * (154.89/4.38277)) #, 'waist' : 355.87 , 'SWEEPS' : 100000 , 'SKIP' : 100 , 'SEED' : seed } )
Write:
h5_infiles = pyalps.writeInputH5Files("parm.nowaist",params);
Or read:
h5_infiles = pyalps.getInputH5Files(prefix = "parm.nowaist");
Parameters (waist)
Parameters:
tof_phase = pyalps.dwa.tofPhase(time_of_flight=15.5, wavelength=[843.,843.,843.], mass=86.99) params=[] for seed in range(10,210,10): params.append( { 'LATTICE' : 'inhomogeneous simple cubic lattice' , 'MODEL' : 'boson Hubbard' , 'L' : 160 , 'Nmax' : 20 , 't' : 1. , 'U' : 8.1 , 'T' : 1. , 'mu_homogeneous' : 2. # Attempt to fix at N=280000 , lower bound: 2 , upper bound: 3 , 'mu' : 'mu_homogeneous - (0.00369*(x-(L-1)/2.)*(x-(L-1)/2.) + 0.00369*(y-(L-1)/2.)*(y-(L-1)/2.) + 0.00369*(z-(L-1)/2.)*(z-(L-1)/2.))' , 'tof_phase' : str(tof_phase) , 'V0' : (8.35 * (154.89/4.38277)) , 'waist' : 355.87 , 'SWEEPS' : 100000 , 'SKIP' : 100 , 'SEED' : seed } )
Write:
h5_infiles = pyalps.writeInputH5Files("parm.waist",params);
Or read:
h5_infiles = pyalps.getInputH5Files(prefix = "parm.waist");
A preliminary run
A trial run like this:
pyalps.runApplication('dwa', h5_infiles[0]);
Automated processing (no waist)
def parallel_automated_processing(taskfile): pyalps.dwa.recursiveRun( "pyalps.runApplication('dwa', taskfile, T=2700)" , cmd_lang = 'python' , follow_up_script = "pyalps.dwa.extract_worldlines(infile=pyalps.input2output(taskfile), outfile=taskfile)" , break_if = "pyalps.dwa.thermalized(pyalps.input2output(taskfile), 'Total Particle Number', tolerance=0.001, simplified=True)" , write_status = "pyalps.dwa.write_status(pyalps.input2output(taskfile), 'Thermalizing')" , loc = locals() , batch_submit = True , batch_cmd_prefix = 'bsub -J jobno -n1 -W60' , batch_run_script = taskfile + '.thermalize.script' , batch_next_run_script = taskfile + '.converge.script' , batch_noRun = True ); pyalps.dwa.recursiveRun( "pyalps.runApplication('dwa', pyalps.input2output(taskfile), T=2700)" , cmd_lang = 'python' , break_if = "pyalps.dwa.converged(pyalps.input2output(taskfile), 'Total Particle Number', simplified=True)" , write_status = "pyalps.dwa.write_status(pyalps.input2output(taskfile), 'Converging')" , loc = locals() , batch_submit = True , batch_cmd_prefix = 'bsub -J jobno -n1 -W60' , batch_run_script = taskfile + '.converge.script' , batch_noRun = True ); pyalps.dwa.startRunScript( batch_run_script= taskfile + '.thermalize.script' , batch_cmd_prefix = 'bsub -J jobno -n1 -W60' , loc = locals() );
for h5_infile in h5_infiles: parallel_automated_processing(h5_infile);
Automated processing (waist)
def parallel_automated_processing(taskfile): pyalps.dwa.recursiveRun( "pyalps.runApplication('dwa', taskfile, T=2700)" , cmd_lang = 'python' , follow_up_script = "pyalps.dwa.extract_worldlines(infile=pyalps.input2output(taskfile), outfile=taskfile)" , break_if = "pyalps.dwa.thermalized(pyalps.input2output(taskfile), 'Total Particle Number', tolerance=0.001, simplified=True)" , write_status = "pyalps.dwa.write_status(pyalps.input2output(taskfile), 'Thermalizing')" , loc = locals() , batch_submit = True , batch_cmd_prefix = 'bsub -J jobyes -n1 -W60' , batch_run_script = taskfile + '.thermalize.script' , batch_next_run_script = taskfile + '.converge.script' , batch_noRun = True ); pyalps.dwa.recursiveRun( "pyalps.runApplication('dwa', pyalps.input2output(taskfile), T=2700)" , cmd_lang = 'python' , break_if = "pyalps.dwa.converged(pyalps.input2output(taskfile), 'Total Particle Number', simplified=True)" , write_status = "pyalps.dwa.write_status(pyalps.input2output(taskfile), 'Converging')" , loc = locals() , batch_submit = True , batch_cmd_prefix = 'bsub -J jobyes -n1 -W60' , batch_run_script = taskfile + '.converge.script' , batch_noRun = True ); pyalps.dwa.startRunScript( batch_run_script= taskfile + '.thermalize.script' , batch_cmd_prefix = 'bsub -J jobyes -n1 -W60' , loc = locals() );
for h5_infile in h5_infiles: parallel_automated_processing(h5_infile);
Tuning chemical potential
Parameters (no waist)
Headers:
import numpy; import pyalps; import pyalps.dwa;
Parameters:
tof_phase = pyalps.dwa.tofPhase(time_of_flight=15.5, wavelength=[843.,843.,843.], mass=86.99) params=[] for mu in [2.0, 2.1, 2.2, 2.3, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9]: params.append( { 'LATTICE' : 'inhomogeneous simple cubic lattice' , 'MODEL' : 'boson Hubbard' , 'L' : 160 , 'Nmax' : 20 , 't' : 1. , 'U' : 8.1 , 'T' : 1. , 'mu_homogeneous' : mu # Attempt to fix at N=280000 , lower bound: 2 , upper bound: 3 , 'mu' : 'mu_homogeneous - (0.00369*(x-(L-1)/2.)*(x-(L-1)/2.) + 0.00369*(y-(L-1)/2.)*(y-(L-1)/2.) + 0.00369*(z-(L-1)/2.)*(z-(L-1)/2.))' , 'tof_phase' : str(tof_phase) #, 'V0' : (8.35 * (154.89/4.38277)) #, 'waist' : 355.87 , 'SWEEPS' : 100000 , 'SKIP' : 100 } )
Write:
h5_infiles = pyalps.writeInputH5Files("parm.nowaist",params);
Or read:
h5_infiles = pyalps.getInputH5Files(prefix = "parm.nowaist");
Parameters (waist)
Headers:
import numpy; import pyalps; import pyalps.dwa;
Parameters:
tof_phase = pyalps.dwa.tofPhase(time_of_flight=15.5, wavelength=[843.,843.,843.], mass=86.99) params=[] for mu in [2.0, 2.1, 2.2, 2.3, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9]: params.append( { 'LATTICE' : 'inhomogeneous simple cubic lattice' , 'MODEL' : 'boson Hubbard' , 'L' : 160 , 'Nmax' : 20 , 't' : 1. , 'U' : 8.1 , 'T' : 1. , 'mu_homogeneous' : mu # Attempt to fix at N=280000 , lower bound: 2 , upper bound: 3 , 'mu' : 'mu_homogeneous - (0.00369*(x-(L-1)/2.)*(x-(L-1)/2.) + 0.00369*(y-(L-1)/2.)*(y-(L-1)/2.) + 0.00369*(z-(L-1)/2.)*(z-(L-1)/2.))' , 'tof_phase' : str(tof_phase) , 'V0' : (8.35 * (154.89/4.38277)) , 'waist' : 355.87 , 'SWEEPS' : 100000 , 'SKIP' : 100 } )
Write:
h5_infiles = pyalps.writeInputH5Files("parm.waist",params);
Or read:
h5_infiles = pyalps.getInputH5Files(prefix = "parm.waist");
Observing on your own
First the headers again:
import matplotlib; matplotlib.use('Agg');
import pyalps; import pyalps.plot; import pyalps.dwa;
We gather result files:
resultFiles = pyalps.getResultFiles(prefix='parm.nowaist');
or
resultFiles = pyalps.getResultFiles(prefix='parm.waist');
and collect data:
data = pyalps.loadMeasurements(resultFiles, 'Total Particle Number'); N = pyalps.collectXY(data, x='mu_homogeneous', y='Total Particle Number')
then plot it:
fig = matplotlib.pyplot.figure(); pyalps.plot.plot(N); fig.savefig('parm.N.jpg'); fig.show();
or send it to yourself:
pyalps.sendmail( 'pingnang@phys.ethz.ch' # email address of recipients , attachment='parm.N.jpg' , subject='Picture: parm.N.jpg' );
Unexpected norun
Sometimes, due to certain node failure, jobs may not even start. The following is the way to fix this:
import pyalps; import pyalps.dwa;
all_h5_infiles = pyalps.getInputH5Files(prefix='parm.nowaist'); h5_infiles = []; for h5_infile in all_h5_infiles: ar = pyalps.hdf5.h5ar(h5_infile); if not ar.is_group('/simulation/worldlines'): h5_infiles.append(h5_infile);
for h5_infile in h5_infiles: parallel_automated_processing(h5_infile);
Checking thermalization
Gather result files:
resultFiles = pyalps.getResultFiles(prefix='parm.nowaist');
or:
resultFiles = pyalps.getResultFiles(prefix='parm.waist');
Check now:
for resultFile in resultFiles: mu=pyalps.getParameters(resultFile)[0]['mu_homogeneous']; N=pyalps.getMeasurements(resultFile,observable='Total Particle Number')[0]['mean']['value']; thermalized=pyalps.dwa.thermalized(resultFile, 'Total Particle Number', tolerance=0.001, simplified=True); print(resultFile + ' , mu : ' + str(mu) + ' , N : ' + str(N) + ' , thermalized : ' + str(thermalized));
Checking convergence
Gather result files:
resultFiles = pyalps.getResultFiles(prefix='parm.nowaist');
or:
resultFiles = pyalps.getResultFiles(prefix='parm.waist');
Check now:
for resultFile in resultFiles: mu=pyalps.getParameters(resultFile)[0]['mu_homogeneous']; N=pyalps.getMeasurements(resultFile,observable='Total Particle Number')[0]['mean']['value']; converged=pyalps.dwa.converged(resultFile, 'Total Particle Number', simplified=True); print(resultFile + ' , mu : ' + str(mu) + ' , N : ' + str(N) + ' , converged : ' + str(converged));
Rethermalization
Gather the files:
h5_infiles = pyalps.getInputH5Files(prefix='parm.nowaist');
or:
h5_infiles = pyalps.getInputH5Files(prefix='parm.waist');
Copy worldlines over:
for h5_infile in h5_infiles: pyalps.dwa.extract_worldlines(pyalps.input2output(h5_infile), h5_infile); print('Done : ' + h5_infile);
Rethermalize:
for h5_infile in h5_infiles: parallel_automated_processing(h5_infile);
Fine tuning the chemical potential
First the headers again:
import numpy; import pyalps; import pyalps.dwa;
Easily done like this:
resultFiles = pyalps.getResultFiles(prefix='parm.nowaist');
or
resultFiles = pyalps.getResultFiles(prefix='parm.waist');
data = pyalps.loadMeasurements(resultFiles, 'Total Particle Number'); params = pyalps.paramsAtFixedY(data, x='mu_homogeneous', y='Total Particle Number', fixedY=280000);
Prepare new input files:
h5_infiles = pyalps.writeInputH5Files("parm.nowaist.1", params);
or
h5_infiles = pyalps.writeInputH5Files("parm.waist.1", params);
Run them:
for h5_infile in h5_infiles: parallel_automated_processing(h5_infile);
Into the business
Parameters (no waist)
Parameters:
tof_phase = pyalps.dwa.tofPhase(time_of_flight=15.5, wavelength=[843.,843.,843.], mass=86.99) params=[] for seed in range(10,210,10): params.append( { 'LATTICE' : 'inhomogeneous simple cubic lattice' , 'MODEL' : 'boson Hubbard' , 'L' : 160 , 'Nmax' : 20 , 't' : 1. , 'U' : 8.1 , 'T' : 1. , 'mu_homogeneous' : 2.54 # Attempt to fix at N=280000 , 'mu' : 'mu_homogeneous - (0.00369*(x-(L-1)/2.)*(x-(L-1)/2.) + 0.00369*(y-(L-1)/2.)*(y-(L-1)/2.) + 0.00369*(z-(L-1)/2.)*(z-(L-1)/2.))' , 'tof_phase' : str(tof_phase) #, 'V0' : (8.35 * (154.89/4.38277)) #, 'waist' : 355.87 , 'SWEEPS' : 100000 , 'SKIP' : 100 , 'SEED' : seed } )
Write:
h5_infiles = pyalps.writeInputH5Files("parm.nowaist",params);
Or read:
h5_infiles = pyalps.getInputH5Files(prefix = "parm.nowaist");
Parameters (waist)
Parameters:
tof_phase = pyalps.dwa.tofPhase(time_of_flight=15.5, wavelength=[843.,843.,843.], mass=86.99) params=[] for seed in range(10,210,10): params.append( { 'LATTICE' : 'inhomogeneous simple cubic lattice' , 'MODEL' : 'boson Hubbard' , 'L' : 160 , 'Nmax' : 20 , 't' : 1. , 'U' : 8.1 , 'T' : 1. , 'mu_homogeneous' : 2.45 # Attempt to fix at N=280000 , 'mu' : 'mu_homogeneous - (0.00369*(x-(L-1)/2.)*(x-(L-1)/2.) + 0.00369*(y-(L-1)/2.)*(y-(L-1)/2.) + 0.00369*(z-(L-1)/2.)*(z-(L-1)/2.))' , 'tof_phase' : str(tof_phase) , 'V0' : (8.35 * (154.89/4.38277)) , 'waist' : 355.87 , 'SWEEPS' : 100000 , 'SKIP' : 100 , 'SEED' : seed } )
Write:
h5_infiles = pyalps.writeInputH5Files("parm.waist",params);
Or read:
h5_infiles = pyalps.getInputH5Files(prefix = "parm.waist");
Turn on respective measurements
Gather the files:
h5_infiles = pyalps.getInputH5Files(prefix='parm.nowaist')
or:
h5_infiles = pyalps.getInputH5Files(prefix='parm.waist')
Switch on local density measurements:
for h5_infile in h5_infiles: pyalps.dwa.switchParameter(h5_infile, 'MEASURE[Local Density]', 1); pyalps.dwa.switchParameter(pyalps.input2output(h5_infile), 'MEASURE[Local Density]', 1);
Switch on green function measurements:
for h5_infile in h5_infiles: pyalps.dwa.switchParameter(h5_infile, 'MEASURE[Green Function]', 1); pyalps.dwa.switchParameter(pyalps.input2output(h5_infile), 'MEASURE[Green Function]', 1);
Automated processing (no waist)
def parallel_automatic_processing(taskfile): pyalps.dwa.recursiveRun( "pyalps.runApplication('dwa', taskfile, T=2700)" , cmd_lang = 'python' , write_status = "pyalps.dwa.write_status(pyalps.input2output(taskfile), 'Measuring')" , loc = locals() , batch_submit = True , batch_cmd_prefix = 'bsub -J jobno -n1 -W60 -r' , batch_run_script = taskfile + '.measure1.script' , batch_next_run_script = taskfile + '.measure.script' , batch_noRun = True ); pyalps.dwa.recursiveRun( "pyalps.runApplication('dwa', pyalps.input2output(taskfile), T=2700)" , cmd_lang = 'python' , n = 20 , write_status = "pyalps.dwa.write_status(pyalps.input2output(taskfile), 'Measuring')" , loc = locals() , batch_submit = True , batch_cmd_prefix = 'bsub -J jobno -n1 -W60 -r' , batch_run_script = taskfile + '.measure.script' , batch_next_run_script = taskfile + '.sendmail.script' , batch_noRun = True ); pyalps.dwa.recursiveRun( "pyalps.sendmail('pingnang@phys.ethz.ch', message= taskfile + ' has finished measurements for 20 iterations.')" , cmd_lang = 'python' , loc = locals() , batch_submit = True , batch_run_script = taskfile + '.sendmail.script' , batch_noRun = True ); pyalps.dwa.startRunScript( batch_run_script= taskfile + '.measure1.script' , batch_cmd_prefix = 'bsub -J jobno -n1 -W60 -r' , loc = locals() );
for h5_infile in h5_infiles: parallel_automatic_processing(h5_infile);
Automated processing (waist)
def parallel_automatic_processing(taskfile): pyalps.dwa.recursiveRun( "pyalps.runApplication('dwa', taskfile, T=2700)" , cmd_lang = 'python' , write_status = "pyalps.dwa.write_status(pyalps.input2output(taskfile), 'Measuring')" , loc = locals() , batch_submit = True , batch_cmd_prefix = 'bsub -J jobyes -n1 -W60 -r' , batch_run_script = taskfile + '.measure1.script' , batch_next_run_script = taskfile + '.measure.script' , batch_noRun = True ); pyalps.dwa.recursiveRun( "pyalps.runApplication('dwa', pyalps.input2output(taskfile), T=2700)" , cmd_lang = 'python' , n = 20 , write_status = "pyalps.dwa.write_status(pyalps.input2output(taskfile), 'Measuring')" , loc = locals() , batch_submit = True , batch_cmd_prefix = 'bsub -J jobyes -n1 -W60 -r' , batch_run_script = taskfile + '.measure.script' , batch_next_run_script = taskfile + '.sendmail.script' , batch_noRun = True ); pyalps.dwa.recursiveRun( "pyalps.sendmail('pingnang@phys.ethz.ch', message= taskfile + ' has finished measurements for 20 iterations.')" , cmd_lang = 'python' , loc = locals() , batch_submit = True , batch_run_script = taskfile + '.sendmail.script' , batch_noRun = True ); pyalps.dwa.startRunScript( batch_run_script= taskfile + '.measure1.script' , batch_cmd_prefix = 'bsub -J jobyes -n1 -W60 -r' , loc = locals() );
for h5_infile in h5_infiles: parallel_automatic_processing(h5_infile);
Local density
The headers:
import matplotlib; matplotlib.use('Agg');
import matplotlib.pyplot; from mpl_toolkits.mplot3d.axes3d import Axes3D; import numpy; import pyalps; import pyalps.dwa;
h5_infiles = pyalps.getInputH5Files(prefix='parm.nowaist'); h5_infile = h5_infiles[0]
L = int(pyalps.getParameters(pyalps.input2output(h5_infile))[0]['L']); density = pyalps.getMeasurements(pyalps.input2output(h5_infile), observable='Local Density')[0]['mean']['value'].reshape([L,L,L]); cross_section_density = density[:,:,L/2]; column_integrated_density = numpy.sum(density, axis=2); x = numpy.array([range(L)] * L, float) - (L-1)/2.; y = numpy.transpose(x);
To plot the cross section density along a line:
matplotlib.pyplot.figure(); matplotlib.pyplot.plot(x[L/2,:], cross_section_density[L/2,:]); matplotlib.pyplot.show();
or column integrated density along a line:
matplotlib.pyplot.figure(); matplotlib.pyplot.plot(y[:,L/2], column_integrated_density[:,L/2]); matplotlib.pyplot.show();
To see a cross section density on a plane:
fig = matplotlib.pyplot.figure(); ax = fig.add_subplot(1, 1, 1, projection='3d') surf = ax.plot_surface( x, y, cross_section_density, rstride=1, cstride=1, cmap=matplotlib.cm.coolwarm, linewidth=0, antialiased=False) fig.colorbar(surf, shrink=0.5, aspect=10); fig.savefig(h5_infile + '.cross-section-density.jpg'); fig.show();
To see a column integrated density on a plane:
fig = matplotlib.pyplot.figure(); ax = fig.add_subplot(1, 1, 1, projection='3d') surf = ax.plot_surface( x, y, column_integrated_density, rstride=1, cstride=1, cmap=matplotlib.cm.coolwarm, linewidth=0, antialiased=False) fig.colorbar(surf, shrink=0.5, aspect=10); fig.savefig(h5_infile + '.column-integrated-density.jpg'); fig.show();
Send the plot to yourself:
pyalps.sendmail( 'pingnang@phys.ethz.ch' # email address of recipients , attachment=h5_infile +'.cross-section-density.jpg' , subject='Picture: ' + h5_infile + '.cross-section-density.jpg' );
pyalps.sendmail( 'pingnang@phys.ethz.ch' # email address of recipients , attachment=h5_infile + '.column-integrated-density.jpg' , subject='Picture: ' + h5_infile + '.column-integrated-density.jpg' );
For this project
The headers:
import matplotlib; matplotlib.use('Agg');
import matplotlib.pyplot; from mpl_toolkits.mplot3d.axes3d import Axes3D; import math; import numpy; import pyalps; import pyalps.dwa;
This:
h5_infiles = pyalps.getInputH5Files(prefix='parm.nowaist');
or this:
h5_infiles = pyalps.getInputH5Files(prefix='parm.waist');
obs = []; for h5_infile in h5_infiles: L = int(pyalps.getParameters(pyalps.input2output(h5_infile))[0]['L']); density = pyalps.getMeasurements(pyalps.input2output(h5_infile), observable='Local Density')[0]['mean']['value'].reshape([L,L,L]); cross_section_1d_density = (density[:,L/2,L/2] + density[L/2,:,L/2] + density[L/2,L/2,:])/3.; cross_section_1d_density = (cross_section_1d_density + cross_section_1d_density[::-1])/2.; x = numpy.array(range(L), float) - (L-1)/2.; obs.append({ 'x' : x , 'obs' : cross_section_1d_density}); print('Done : ' + h5_infile);
x = obs[0]['x']; obs = numpy.array([item['obs'] for item in obs]); value = obs.mean(axis=0); error = obs.std(axis=0)/math.sqrt(obs.shape[0]);
Data
Saving:
numpy.savetxt('parm.nowaist.cross-section-1d-density.dat', numpy.array([x,value,error]).transpose())
or:
numpy.savetxt('parm.waist.cross-section-1d-density.dat', numpy.array([x,value,error]).transpose())
Email:
Send the plot to yourself:
pyalps.sendmail( 'pingnang@phys.ethz.ch' # email address of recipients , attachment='parm.nowaist.cross-section-1d-density.dat' , subject='Picture: parm.nowaist.cross-section-1d-density.dat' );
or
pyalps.sendmail( 'pingnang@phys.ethz.ch' # email address of recipients , attachment='parm.waist.cross-section-1d-density.dat' , subject='Picture: parm.nowaist.cross-section-1d-density.dat' );
Plotting figure
To plot the cross section density along a line:
fig = matplotlib.pyplot.figure(); matplotlib.pyplot.errorbar(x, value, error);
fig.savefig('parm.nowaist.cross-section-1d-density.jpg');
or
fig.savefig('parm.waist.cross-section-1d-density.jpg');
Send the plot to yourself:
pyalps.sendmail( 'pingnang@phys.ethz.ch' # email address of recipients , attachment='parm.nowaist.cross-section-1d-density.jpg' , subject='Picture: parm.nowaist.cross-section-1d-density.jpg' );
or
pyalps.sendmail( 'pingnang@phys.ethz.ch' # email address of recipients , attachment='parm.waist.cross-section-1d-density.jpg' , subject='Picture: parm.nowaist.cross-section-1d-density.jpg' );