
MMJCO - Tunnel Junction Amplitude calculator.
Version 0.95 (ALPHA+)    08/20/2022
S. Whiteley,  wrcad.com, Synopsys, Inc.


WRspice contains an internal tunnel junction model (TJM) of a
Josephson junction.  The model requires pre-prepared tables of "fit"
parameters, which avoid performing lengthly calculation at run time. 
The mmjco program generates the needed files.

In general, the user may not need to interact directly with mmjco, as
WRspice will call mmjco when needed to create new fit files, and store
these under a directory named ".mmjco" in the user's home directory. 
Once a fit file has been created for a particular parameter set, it
will be reused in future WRspice sessions when needed,

The core functionality of mmjco is derived from the MiTMoJCo project
by D.  R.  Gulevich (found on GitHub.com), specifically the
mitmojco.py environment that provides an interface for creating tunnel
junction amplitude (TCA) and fitting tables.  This has been
implemented in mmjco> as two C++ classes, one for creating tunnel
junction amplitudes, which basically evaluates the Werthamer model as
formulated for MiTMoJCo by Gulevich, the second to compress the
amplitudes into a compact representation.  This representation is used
by the TJM (JJ level=3) device model in WRspice.  The method is that
of Odintsov, Semenov, and Zorin (See below for a list of references).

The <b>mmjco</b> program has additional features.
o  Sweep tables can be produced for a range of temperatures, allowing
   temperature sweeps to be performed in simulation without having to
   create a fit file at each temperature.

o  There is a built-in BCS computation of energy gap given
   tempoerature, superconducting transition temperature, and Debye
   temperature of the junction materials.

o  Tunnel current amplitude tables are generated in "rawfile" format,
   so can be plotted within WRspice or Synopsys WaveView.

o  The fit and sweep file formats are compatible with the TJM developed
   for Synopsys PrfimeSym-HSPICE under the IARPA SuperTools project.


++++ Building mmjco

One should be able to build mmjco on any Linux-like platform.

Presently, the file will build into a static library, and a
stand-alone application "mmjco" that will provide an interactive
interface for generating files containing tunnel junction amplitudes
and fitting tables.

Building requires the "devel" versions of the following libraries:

    GSL  (the Gnu Scientific library)  Version 1.15 is known to be ok.

The source now includes a release of cmpfit
https://pages.physics.wisc.edu/~craigm/idl/cmpfit.html
This provides a version of the venerable lmdif least-squares fitter
from MINPACK, an ancient fortran library.  As an alternative, one can
instead use:
    cminpack (a C/C++ version of the MINPACK optimization library)
    https://github.com/devernay/cminpack

This was the original target for mmjco, however it is a bit obscure,
and in particular I couldn't find it at Synopsys (though it is
provided in CentOS 7 repos).  This may be an issue for others so I
found and included cmpfit instead.

Once installed, tweak the Makefile if necessary to set the proper
paths to libraries for your system, and type "make".  If all goes well
the "mmjco" binary should be created in the current directory, along
with a libmmjco.a library file.

If built with the preprocessor symbol WRSPICE defined, the program is
intended to serve as an accessory to the WRspice circuit simulator,
for building tables for use with the TMJ model.  There are two main
effects:

  - The tunnel current amplitude files use the "rawfile" format, and
    can be plotted in WRspice, WaveView, and other programs.

  - The created amplitude and fit files are saved by default in a
    directory named ".mmjco" in the user's home directory.  This
    location is in the default WRspice search path for fit files.


++++ Running mmjco

Running mmjco enters a command-line processing loop, the user responds
to the "mmjco> " prompt with commands from the list below.  Each
command can be followed by options as indicated.  Additionally,there
are "cdf" and "swp" modes where mmjco will create TCA and fit files or
temperature sweep files according to the arguments, and exit.  This is
mostly to support WRspice, which uses this mode to create new models
on-the-fly.

In this document, text in square brackets ([...]) is optional.  A
vertical pipe character ( | ) indicates that either the text to the
right or left is acceptable, i.e., it separates options.


++++ Command Line Operations

mmjco cdf arguments

If the first argument to the mmjco executable is "cdf", a TCA file and
corresponding fit file are created, and mmjco exits immediately in
this case.  Arguments following "cdf" are the same arguments that can
follow the cd and cf interactive commands.  This mode is used by the
TJM device model in WRspice.

mmjco swp -fs sweepfile temperature

Similarly, "swf" will create a possibly-interpolated fit file
from an existing sweep file, for the temperature provided.


++++ Interactive mmjco Commands

The following are the interactive commands which can be entered after
starting mmjco with no arguments.

cd[ata]  [-t temp] [-d|-d1|-d2 delta] [-s smooth] [-x nx] [-f filename] \
         [-r | -rr | -rd]

    Create TCA data, save internally and to a file.  See below for
    an explanation of the options.

    Tunnel current amplitude and smoothing options:

    -t   The assumed temperature follows, in Kelvin.  Default 4.2.

    -d   This will set both d1 and d2, the pair breaking energy in
         milli-electron volts, of the two superconducting banks.  The
         default is 1.4 mev.

    -d1  
    -d2  Like -d, but apply to only one of the banks.  The final
         occurrence of d,d1,d2 will have precedence.

    -s   This provides the smoothing parameter as used in MiTMoJCo.
         the accepted range is 0.0 - 0.099.  The default is 0.008.  If
         less than 0.001, 0 is assumed.  When 0, no smoothing is done
         and raw BCS tunnel amplitudes are generated.

    -x   The number of points used to create the tunnel current
         amplitudes.  The range of sweep of voltage normalized to the
         gap voltage (d1+d2) extends from 0.001 through 2.0.  The
         default point count is 500.

    -f   A name for the TCA anplitude file.  If not given, a
         default is used, described below.

    -r   Output file is a complex-valued rawfile.
    -rr  Output file is a real-valued rawfile.
    -rd  Output file diomple data file.  If none of the "-r" options
         is set, the format will be "-rr" if the program was built with
         WRSPICE defined, "-rd" otherwise.

cf[it]  [-n terms] [-h thr] [-ff filename]

    Create fit parameters for TCA data currently in memory, from cd
    or ld commands.  This is saved internallyand to a file.  See
    below for an explanation of the options.

    Fitting table options:

    -n   The size of the table, defaults to 8.  Larger tables are
         more accurate but take more time to generate and process.
         A maximum of 20 is enforced.

    -h   The ratio of the absolute to relative tolerances, used in
         compression, the default is 0.2.

    -ff  A name for the fitting parameter table.  If not given, a
         default is used, described below.

cm[odel]  [-h thr] [-fm [filename]] [-r | -rr | -rd]

    Create a model for TCA data using fitting parameters currently in
    memory, compute the residual, and optionally save to a file.  If -fm
    is given without a filename, a file name will be generated
    internally.  If -fm is not given, the model will not be saved to a
    file, but used only to compute the residual.  The printed
    residual number is an indication of the fit quality, smaller
    values indicate better matching.

    If one of -r, -rr, -rd is given when a file is being generated, it
    overrides the default type for file to produce.
      -r  complex-valued rawfile
      -rr real-valued rawfile
      -rd simple data file

    The default file type is "-rr" when built with WRSPICE defined,
    or the simple data file format otherwise.

    The model files are saved in the current directory (unless a path
    is given explicitly).


cs[weep] Tstrt Tend [Tdelta] arguments

    Create a temperature sweep file, which involves creating
    sequential records of fit parameters for temperatures starting
    with Tstrt and ending at at or near Tend, spaced in temperature by
    Tdelta.  These are real numbers in Kelvin.  If the third numbver
    Tdelta does not appear, 0.1K is assumed.  The arguments are those
    that can be given to the cd or cf commands.

ct[ab] T1 T2 [... TN] arguments

    Create a temperature table file.  The first two temperatures T1
    and T2 are required, and these can be followed by an arbitrary
    number of additional temperatures.  The temperatures are real
    numbers in Kelvin.  The arguments are those that can be given to
    the cd or cf commands.  The file will contain fit parameters for
    each temperature given.

d[ir] directory_path

    Give a path to a directory where all amplitude and fit files will
    be stored and loaded from, if a rooted path is not given with the
    file names.  The location is by default the current directory,
    unless the WRSPICE preprocessor symbol is used, in which case the
    default location is ".mmjco" in the user's home directory.

g[ap] [-tc Tc] [-td Td] [T1 T2...]

    Compute and print the superconducting energy gap at temperatures. 
    The -tc specifies the superconducting transition temperature, and
    -td specifies the Debye temperature, both Kelvin.  If not given,
    the defaults are for Niobium:  Tc = 9.26K and Td = 276K.  There
    can be zero to 10 real number arguments representing temperatures
    in Kelvin.  If zero, print the gap for temperatures from 0 to Tc
    at 0.1K increments.  If two numbers, print the gap for the smaller
    to the larger in 0.1K increments.  Otherwise, print the gap at the
    given temperatures.

ld[ata] filename

    Load the internal TCA data register from a TCA data file whose name
    must be given.  This understands all supported file formats.

lf[it] filename

    Load the internal fit parameters register from a fit parameter file
    given.

ls[weep] -fs filename temp

    Load the internal fit register from a sweep file, interpolating to
    temperature temp.

h[elp] | v[ersion] | ?

    Print help and the running mmjco release number.

q[uit] | e[xit]

    Exit mmjco.


++++ File Name Encoding

Running mmjco can generate three types of files:  a file containing
values of the tunnel current amplitudes (TCA files), a file containing
the fitting parameters used to represent the tunnel current amplitues
in a compact form, and files that contain collections of fitting
parameter sets, as for a temperature sweep.  The "model" files use the
same format as the TCA files.  By default, these files are given a
name that encodes the various parameter values used in creation.

Tunnel current amplitude files, including model files, are by default
given an internally-generated file name in the forms below.

    tcaTTTTTTdddddDDDDDssPPPP.data
    tcaTTTTTTdddddDDDDDssPPPP.raw

The "tca" and ".data"/".raw" and similar are literal.  All fields are
fixed width and zero padded.  Real numbers are converted to integers
by multiplying by a scale factor and rounding to integer.

TTTTTT
    A six-digit integer equal to 1e4*<i>temperature</i>.
ddddd
    A 5-digit integer equal to 1e7*<i>d1</i>, where <i>d1</i> is the
    side 1 gap potential in volts.
DDDDD
    A 5-digit integer equal to 1e7*<i>d2</i>, where <i>d2</i> is the
    side 2 gap potential in volts.
ss
    A two digit integer equal to 1e3*<i>sf</i>, where <i>sf</i> is the
    smoothing parameter.
PPPP
    The number of points used (<tt>-x</tt> option).

In WRspice, a rawfile (.raw extension) can be loaded with the load
command, and the amplitudes can then be plotted with "plot all".

The fitting parameter file name has the following form.

    tcaTTTTTTdddddDDDDDssPPPP-NNHHH.fit

The first part, ahead of the '-', is as described above.  Following
the hyphen:

NN
    The two-digit fitting table size.
HHH
    A three digit integer equal to 1e3*thr, where thr is the
    compression threshold.

Sweep file names have the following form.

    tswNNNttttttTTTTTTssPPPP-NNHHH.swp

The part following the hyphen is the same as for the fit file.

NNN
    A three digit integer giving the number of records in the file.
tttttt
    A six digit integer equal to 1e4*Tstart, where Tstart is the
    starting temperature in Kelvin.
TTTTTT
    A six digit integer equal to 1e4*Tdelta, where Tdelta is the
    temperature spacing.
ss
    A two digit integer equal to 1e3*sf, where sf is the
    smoothing parameter.
PPPP
    The number of points used (-x option).

The fit parameter table format is similar to that used by MiTMoJCo,
identical if the header is ignored.  Note that the table values do not
match those found in the MiTMoJCo distribution.  It seems that these
values are not unique, and the various programs can converge to
different sets.  It was found that the original and present versions
of mitmojco.py gave different values, and neither matched the values
provided in the amplitudes folder.


++++ File Formats

TCA file formats
    Tunnel current amplitude (TCA) tables are created by the cd and cm
    commands.  They consist of real and imaginary parts of the pair
    and quasiparticle amplitudes, on a scale of normalized potential
    across the structure, where the unit value is the sum of the gap
    energies of the two electrodes.  The scale extends from 0 to 2 in
    these units.

    There actually are three formats available, selectable with
    options to the cd and cm commands.

    option: -rd
    The file name is described in the previous section, with a suffix
    ".data".  This is a generic numerical format, consisting of a
    comment line and six columns of numbers.  The number of data lines
    is the value given with the "-x" option to the cd command,
    defaulting to 500.  The example below illustrates the format.

#    X            Jpair_real   Jpair_imag   Jqp_real     Jqp_imag
0    1.00000e-03  7.50623e-01  2.85275e-04  -7.50625e-01 3.37356e-04
1    5.00601e-03  7.51136e-01  1.36358e-03  -7.51123e-01 1.62410e-03
2    9.01202e-03  7.51643e-01  2.27465e-03  -7.51594e-01 2.74301e-03
3    1.30180e-02  7.52144e-01  3.03372e-03  -7.52040e-01 3.70900e-03
...
498  1.99599e+00  4.22402e-01  -5.29166e-01 1.52235e-02  1.87539e+00
499  2.00000e+00  4.21420e-01  -5.28564e-01 1.51180e-02  1.87962e+00

    The other two options emit the data using the SPICE "rawfile"
    format.  This is a format developed for plot data in Berkeley
    Spice3, which is supported by most plotting programs, including
    Synopsys WaveView and the load function of WRspice.  The only
    difference is that one format ouptuts complex numbers for two
    variables (pair and quasiparticle amplitudes), while the other
    format outputs real values for four variables (the real and
    imaginary parts of the amplitudes).

    option: -r
    The file name is described in the previous section, with a suffix
    ".raw".  Output is in rawfile format using complex numbers.

    option: -rr
    The file name is described in the previous section, with a suffix
    ".raw".  Output is in rawfile format using real numbers.

    The rawfile format description can be found in the WRspice
    documentation.

Fit file format
    A fit file contains a compacted digest of a TCA table, as
    prescribed by the Odintsov, Semenov and Zorin (OSZ) algorithm. 
    These can be generated with the cf command.  Fit files can be used
    as input to simulators that contain a compatible tunnel junction
    model (TJM).  Presently, WRspice and Synopsys HSPICE can use these
    files.

    An example fit file is shown below.

tcafit  4.2000e+00  1.3696e-03  1.3696e-03 0.008 500  8  0.200 3.9580e-2
-5.55136e+00, 7.35249e-02, 1.28573e+00, 1.08362e+01,-1.58776e-01, 2.87279e+01
-1.24623e-02, 1.00032e+00, 5.84894e-03,-4.50828e-04, 5.81876e-03,-2.11396e-04
-3.83312e-02, 1.00112e+00, 2.10779e-02,-3.09528e-04, 2.13882e-02, 1.13197e-03
-1.18261e-01, 9.98252e-01, 6.29481e-02, 1.95061e-03, 6.04875e-02, 1.68472e-02
-5.41049e-02, 7.52572e-07,-2.72420e-04, 5.67468e+00, 5.43002e-04, 2.37802e+00
-9.94491e-01, 6.50936e-01, 7.80341e-01, 1.14849e-01,-2.36003e-01, 9.49693e-01
-3.42835e-01, 9.60836e-01, 1.75970e-01, 1.81214e-02, 1.47672e-01, 1.40112e-01
-2.80018e-01, 6.49039e-03, 6.26863e-03, 1.81423e-01, 9.20013e-03, 3.23103e-02

    The first line is a header, the first word of which is "tcafit".
    The numbers that follow in this line are:
      - The temperature in Kelvin (4.2000e+00).
      - The left electrode pair-breaking energy in ev (1.3696e-03).
      - The right electrode pair-breaking energy in ev (1.3696e-03).
      - The smoothing parameter value used to create the TCA
        table, -s option in the cd command for example (0.008).
      - The number of scale points used in the TCA table, -x option
        in the cd command for example (500).
      - The number of terms used in the fit table, -n option in the
        cf command (8).
      - The value of the threshold parameter used when generating the
        fit parameter table, -h option of the cf command (0.200).
      - Normalized quasiparticle current at x=0.8, used to estimate the
        sub-gap conductance (3.9580e-2).

    Following the header line are six columns of real numbers.  The
    number of rows is equal to the "terms", which is the -n option to
    the cf command.  The columns are the OSZ parameters P.real, P.imag,
    A.real, A.imag, B.real, B.imag.

Sweep file format
    Temperature sweep files are concatenations of fit records as
    described above for a temperature range.  These allow rapid
    temperature modeling through interpolation in supporting
    simulators (WRspice and Synopsys HSPICE).  Sweep files are created
    with the cs command.

    Below is an example temperature sweep file.

tsweep 91 0.1000 0.1000 0.008 500 8 0.200
tcafit  1.0000e-01  1.4086e-03  1.4086e-03  1.3185e-02
-8.51331e+00, 1.15164e-01, 1.29900e+00, 1.11105e+01,-4.05570e-01, 5.65607e+01
-1.06700e-02, 1.00009e+00, 3.58651e-03,-2.25268e-04, 3.57013e-03,-5.66586e-05
-2.47980e-02, 1.00072e+00, 1.14692e-02,-5.49747e-04, 1.15840e-02,-1.07382e-04
-6.31470e-02, 1.00196e+00, 2.92347e-02,-1.70767e-03, 2.93324e-02, 2.09767e-03
-1.86179e+00, 9.23541e-01, 7.52628e-01,-8.71131e-02,-4.37397e-01, 2.11473e+00
-1.56178e-01, 1.00378e+00, 6.81933e-02,-6.76474e-03, 7.04144e-02, 1.46763e-02
-3.74403e-01, 1.01258e+00, 1.53951e-01,-3.71976e-02, 1.73747e-01, 7.86682e-02
-8.73847e-01, 1.06711e+00, 2.56553e-01,-1.53780e-01, 4.86150e-01, 3.82557e-01
tcafit  2.0000e-01  1.4086e-03  1.4086e-03  1.3185e-02
-8.51331e+00, 1.15164e-01, 1.29900e+00, 1.11105e+01,-4.05570e-01, 5.65607e+01
-1.06700e-02, 1.00009e+00, 3.58651e-03,-2.25268e-04, 3.57013e-03,-5.66586e-05
---

    The first line is a file header starting with the word "tsweep". 
    The numbers that follow on this line are:
      - The number of fit records contained in this file (91).
      - The lowest temperature K used for fit parameters in the file,
        this will be used in the first fit record (0.1000).
      - The temperature delta K used in the sweep file (0.1000).
      - The smoothing parameter, -s option, used to create all TCA
        tables (0.008).
      - The number of scale points, -x option, used to create all TCA
        tables (500).
      - The number of terms used for each fit table, -n option (8).
      - The value of the threshold parameter used in each fit table (0.200).

    Following this header, fit records are concatenated.  These are
    similar to the format described above, the only difference is that
    the header line is simplified to omit redundant information.  The
    fit record header contains the following value following the word
    "tcafit".
      - The temperature in Kelvin (4.2000e+00).
      - The left electrode pair-breaking energy in ev (1.3696e-03).
      - The right electrode pair-breaking energy in ev (1.3696e-03).
      - Normalized quasiparticle current at x=0.8, used to estimate the
        sub-gap conductance (3.9580e-2).


++++ References

Tunnel current calculation:
A. I. Larkin and Yu. N. Ovchinnikov, Sov. Phys. JETP 24, 1035 (1967).
D. R. Gulevich, V. P. Koshelets, and F. V. Kusmartsev, Phys. Rev. B 96,
 024515 (2017).
A. B. Zorin, I. O. Kulik, K. K. Likharev, and J. R. Schrieffer, Sov. J. 
 Low Temp. Phys. 5, 537 (1979).
D. R. Gulevich, V. P. Koshelets, and F. V. Kusmartsev, Phys.  Rev. B 96,
 024515 (2017).
D. R. Gulevich, V. P. Koshelets, F. V. Kusmartsev, arXiv:1709.04052 (2017).
D. R. Gulevich, L. V. Filippenko, V. P. Koshelets, arXiv:1809.01642 (2018).

Compression:
A. A. Odintsov, V. K. Semenov and A. B. Zorin, IEEE Trans. Magn. 23, 763
 (1987).
D. R. Gulevich, V. P. Koshelets, and F. V. Kusmartsev, Phys. Rev. B 96,
  024515 (2017).

