gpecn1Dtdpbc.sh

#!/bin/bash
# ======================================================================
# NAME
#
#     gpecn1Dtdpbc.sh
#
# DESCRIPTION
#
#     A bash script to run gpecn1Dtdpbc.x and its supporting programs to
#     complete a full simulation.
#
# NOTES
#
#     All dimensional quantities defined here are assumed to be in terms 
#     of simple harmonic oscillator units.
#
# AUTHOR
#
#     Marty Kandes, Ph.D.
#     Computational & Data Science Research Specialist
#     High-Perfomance Computing User Services Group
#     San Diego Supercomputer Center
#     University of California, San Diego
#
# COPYRIGHT
#     
#     Copyright (c) 2010 - 2021 Martin Charles Kandes
#
# LAST UPDATED
#
#     Tuesday, April 13th, 2021
#
# ----------------------------------------------------------------------

# Declare the version number of gpecn1Dtdpbc.
declare -r version_number='0.5.8'

# Declare system name the simulation will be performed on.
declare -r system_name='expanse'

# Declare run type.
declare -r run_type='sagnac'

# Declare run date.
declare -r run_date='20210410'

# Declare run number.
declare -r run_number='001'

# Define run directory path.
declare -r run_directory="gpecn1Dtdpbc-v${version_number}-${system_name}-${run_type}-${run_date}-${run_number}"

# Define total number grid points on the one-dimensional ring (S1).
declare -r number_of_grid_points=6283

# Define total number of time steps in simulation.
declare -r number_of_time_steps=2000000

# Define the number of time steps to skip between writing output to disk.
declare -r number_of_time_steps_between_writes=2000

# Define the number of Crank-Nicolson iterations performed per time step.
declare -r number_of_iterations_per_time_step=1

# Always define the gamma_factor >= 1.0. This parameter sets the CFL 
# condition for the simulation to maintain stability and prevent 
# suprious oscillations from arising in the solutions.
declare -r gamma_factor=1.0

# Define initial state wave function. 
#
#   wave_function_switch=0
#
#     wavefunc.x generates a (Gaussian) ground state solution of a 
#     one-dimensional simple harmonic oscillator (SHO) potential and 
#     places it into a coherent superposition of two momentum states to 
#     create two propagating Gaussian wave packets.
#
#   wave_function_switch=1
#
#     wavefunc.x generates a superposition of two quantized plane wave 
#     solutions for a one-dimensional ring.
#
#   wave_function_switch=2
#
#     wavefunc.x generates a pair of Bogoliubov-like plane wave
#     excitations on top of a background ground-state solution for a
#     one-dimensional ring.
#
declare -r wave_function_switch=0

# Define external potential acting on wave function. 
#
#   external_potential_switch=0
#
#     potential.x generates a zero potential for all points on the ring.
#
#   external_potential_switch=1
#
#     potential.x generates a one-dimensional simple harmonic oscillator
#     potential centered around an angular position on the ring.
#
declare -r external_potential_switch=0

# Define nonlinear coupling.
#
#   nonlinear_coupling_switch=0
#
#     nonlinear.x generates a uniform, time-independent, and real-valued
#     coupling constant that parameterizes the strength of a (simple, 
#     single-component) Bose-Einstein condensate’s nonlinear mean-field
#     interaction arising from two-body contact-like interactions 
#     between atoms within the condensate.
#
#     nonlinear.x was originally intented to explore many different
#     spatially-complex nonlinear couplings, some even with imaginary
#     components. However, this work was not pursused. As such, no other 
#     options for nonlinear_coupling_switch are currently available.
#
declare -r nonlinear_coupling_switch=0

# Set real/imaginary time mode.
#
#   imaginary_time_switch=0
#
#     Always set imaginary_time_switch=0. gpecn1Dtdpbc.x was originally
#     intended to support impaginary time propagating (ITP) to generate
#     initial state wave function solutions. However, this feature was 
#     never implemented. Only analytic initial state wave functions may
#     be used.
#
declare -r imaginary_time_switch=0

# Define the radius of the ring
declare -r radius_of_ring=10.0

# Define oscillator strength of Gaussian wave packets and/or simple 
# harmonic oscillator potential.
declare -r gaussian_oscillator_strength=0.1

# Define initial angular position of initial wave packets generated by 
# wave_function_switch=0 and/or simple harmonic oscillator potential 
# created by external_potential_switch=1.
declare -r initial_angular_position=1.570796326794897

# Define amplitude of first component in the two momentum state 
# superpositions created by wavefunc.x's wave_function_switch=0,1 
# options. The amplitude of the second component is simply defined as 
# the complement of amplitude_1 to conserve (unit) probability. i.e., 
# amplitude_2=sqrt(1-amplitude_1^2). Therefore, for a symmetric 
# superposition set amplitude_1 = 1/sqrt(2).
declare -r amplitude_1=0.707106781186548
#declare -r amplitude_1=0.070710678

# When the wave_function_switch=2, both amplitude_1 and amplitude_2 are
# set independently. Most importantly, however, the sum of the squares 
# of both values must be less than 1. In general, it is assumed that 
# these values will be much less than 1 in order to generate 
# small-amplitude Bogoliubov-like plane wave excitations on a constant
# background ground-state. For a symmetric plane wave superposition on
# top of the ground-state, set amplitude_1 = amplitude_2.
declare -r amplitude_2=0.707106781186548
#declare -r amplitude_2=0.070710678

# Define initial phases associated with the two quantized plane wave 
# solutions in the superposition created by wave_function_switch=1,2.
declare -r initial_phase_1=0.0
declare -r initial_phase_2=0.0

# Define initial angular momenta imparted to the two momentum states 
# created in wavefunc.x's wave_function_switch=0,1,2 options.
declare -r angular_momentum_1=10
declare -r angular_momentum_2=-10

# Define carrier frequency analyzed by ftpsd's Fourier transform phase 
# shift determination algorithm. See K. A. Goldberg and J. Bokor. 
# Fourier-Transform Method of Phase-Shift Determination. Applied Optics,
# 40(17):2886–2894, 2001. 
declare -r carrier_frequency="$(( ${angular_momentum_1} - ${angular_momentum_2} ))"

# Define strength of nonlinear mean-field interaction.
declare -r nonlinear_coupling_amplitude=0.0

# Define angular velocity of the rotating one-dimensional ring.
declare -r rotation_rate_of_ring=0.001

# Make run directory.
mkdir -p "${PWD}/${run_directory}"
cp "${PWD}/wavefunc.x" "${PWD}/${run_directory}"
cp "${PWD}/potential.x" "${PWD}/${run_directory}"
cp "${PWD}/nonlinear.x" "${PWD}/${run_directory}"
cp "${PWD}/gpecn1Dtdpbc.x" "${PWD}/${run_directory}"
cp "${PWD}/expectation.x" "${PWD}/${run_directory}"
cp "${PWD}/ftpsd.x" "${PWD}/${run_directory}"
cp "${PWD}/binary.x" "${PWD}/${run_directory}"
cp "${PWD}/momentum.x" "${PWD}/${run_directory}"
cp "${PWD}/spacetime.x" "${PWD}/${run_directory}"
cp "${PWD}/freqtime.x" "${PWD}/${run_directory}"
cd "${PWD}/${run_directory}"

echo "$(set -o posix; set)" > "${run_directory}.output"

# Compute initial state wave function.
./wavefunc.x \
  "${number_of_grid_points}" \
  "${wave_function_switch}" \
  "${radius_of_ring}" \
  "${rotation_rate_of_ring}" \
  "${gaussian_oscillator_strength}" \
  "${initial_angular_position}" \
  "${amplitude_1}" \
  "${amplitude_2}" \
  "${angular_momentum_1}" \
  "${angular_momentum_2}" \
  "${initial_phase_1}" \
  "${initial_phase_2}" > wavefunc.output

# Compute time-independent external potential.
./potential.x \
  "${number_of_grid_points}" \
  "${external_potential_switch}" \
  "${radius_of_ring}" \
  "${gaussian_oscillator_strength}" \
  "${initial_angular_position}" > potential.output

# Compute nonlinear, mean-field potential.
./nonlinear.x \
  "${number_of_grid_points}" \
  "${nonlinear_coupling_switch}" \
  "${radius_of_ring}" \
  "${nonlinear_coupling_amplitude}" > nonlinear.output

# Simulate.
./gpecn1Dtdpbc.x \
  "${number_of_grid_points}" \
  "${number_of_time_steps}" \
  "${number_of_time_steps_between_writes}" \
  "${number_of_iterations_per_time_step}" \
  "${imaginary_time_switch}" \
  "${radius_of_ring}" \
  "${rotation_rate_of_ring}" \
  "${gamma_factor}" > gpecn1Dtdpbc.output

# Compute expectation values.
./expectation.x \
  "${number_of_grid_points}" \
  "${number_of_time_steps}" \
  "${number_of_time_steps_between_writes}" \
  "${radius_of_ring}" \
  "${rotation_rate_of_ring}" \
  "${gamma_factor}" > expectation.output

# Write (probability) density, real, and imaginary components of binary
# wave function to ascii formatted file
start_file=1000
end_file="$(( ${start_file}+${number_of_time_steps}/${number_of_time_steps_between_writes} ))" 
current_file="${start_file}"
output_shift=10000
while [[ "${current_file}" -le "${end_file}" ]]; do
  ./binary.x \
    "${number_of_grid_points}" \
    "${current_file}" \
    "${output_shift}" \
    "${radius_of_ring}" >> binary.output
  ((current_file+=1))
done

# Compute phase shifts from wave function interference over a range of 
# frequency/momentum components
sideband_frequency_span=10
start_frequency="$(( carrier_frequency - sideband_frequency_span ))"
end_frequency="$(( carrier_frequency + sideband_frequency_span ))"
for (( current_frequency = start_frequency; \
       current_frequency <= end_frequency; \
       current_frequency++ )); do
  ./ftpsd.x \
    "${number_of_grid_points}" \
    "${number_of_time_steps_between_writes}" \
    "${start_file}" \
    "${end_file}" \
    "${radius_of_ring}" \
    "${gamma_factor}" \
    "${current_frequency}" >> ftpsd.output
   mv fort.500 fort.500."${current_frequency}"
done

# Compute momentum distribution.
current_file="${start_file}"
output_shift=20000
minimum_momentum=-100
maximum_momentum=100
while [ "${current_file}" -le "${end_file}" ]; do
  ./momentum.x \
    "${number_of_grid_points}" \
    "${current_file}" \
    "${output_shift}" \
    "${minimum_momentum}" \
    "${maximum_momentum}" \
    "${radius_of_ring}" >> momentum.output
  ((current_file+=1))
done

# Create space-time data.
output_file=591
./spacetime.x \
  "${number_of_grid_points}" \
  "${number_of_time_steps}" \
  "${number_of_time_steps_between_writes}" \
  "${start_file}" \
  "${end_file}" \
  "${output_file}" \
  "${radius_of_ring}" \
  "${gamma_factor}" > spacetime.output

# Create freq-time data.
output_file=592
./freqtime.x \
  "${number_of_grid_points}" \
  "${number_of_time_steps}" \
  "${number_of_time_steps_between_writes}" \
  "${start_file}" \
  "${end_file}" \
  "${output_file}" \
  "${minimum_momentum}" \
  "${maximum_momentum}" \
  "${radius_of_ring}" \
  "${gamma_factor}" > freqtime.output

rm wavefunc.x
rm potential.x
rm nonlinear.x
rm gpecn1Dtdpbc.x
rm expectation.x
rm ftpsd.x
rm binary.x
rm momentum.x
rm spacetime.x
rm freqtime.x

cd ../
tar -czf "${run_directory}".tar.gz "${run_directory}"
rm -rf "${run_directory}"

# ======================================================================