# RAPTOR Tutorial 1: Simplest possible RAPTOR run with restructured procedure of RAPTOR model building.

RAPTOR (RApid Plasma Transport simulatOR) is a code that rapidly solves the tokamak 1D plasma transport PDEs.

It has been designed to be modular and simple to run from a single script, by modifying parameters.

## Contents

- RAPTOR workflow
- The simplest possible RAPTOR run

## RAPTOR workflow

A typical stand-alone run of RAPTOR consists of the steps:

1. Define environment (add path)
2. load default configuration
3. Customize default configuration
4. Create RAPTOR model and parameters for custom configuration.
5. Define input variables (actuator time trajectories)
6. Prepare initial conditions
7. Run RAPTOR.

## The simplest possible RAPTOR run

Before going into details, we will perform all the above steps once in their simplest form.

First, add the relevant paths

    close all hidden;
    run(fullfile(pwd,'..','RAPTOR_path.m')); % add RAPTOR path

Then, load the default configuration structure. This structure contains the default settings from which to configure a RAPTOR simulation

    [config] = RAPTOR_config; % load default params

This loads the default configuration structure, representing a TCV-type plasma setup. We could optionally modify this structure here, but for now we will use the default structure directly. Generate the following structures for these params: model, params, init, geometry g, kinetic profiles v.

    [model,params,init,g,v,U] = build_RAPTOR_model(config);

This generates a RAPTOR model, parameter, and initialization structure. as well as a g and v vector, or matrix.

1. The model structure contains everything that can be pre-calculated, including matrices with spline bases, fixed profiles, and much more. This should not be changed manually but only generated from a call of `RAPTOR_config`
2. The params structure contains everything that can be changed between two runs.
3. The init structure contains information used to define the initial condition in the `RAPTOR_initial_conditions` function (see later)
4. g contains information about geometric quantities related to the equilibrium, which may or may vary in time.
5. v contains information about those kinetic profiles and other quantities which are not in the state (not solved for) and are pre-assigned, and also H-mode and sawtooth crash information.
6. U contains the actuator inputs in time.

    % Define initial conditions for state (depends for example on Ip0 (u(1,1)
    % and other inputs
    % * x is the state vector, which is updated by solving the PDEs.
    % So should be called last before RAPTOR_predictive
    x0 = RAPTOR_initial_conditions(model,init,g(:,1),v(:,1));

Now we are ready to run RAPTOR - this will display some plasma quantities as they evolve.

    simres = RAPTOR_predictive(x0,g,v,U,model,params);

    % To compute some physical outputs, we run the RAPTOR_out function.
    out = RAPTOR_out(simres,model,params);

         it  telaps newt     res   t[ms]  dt[ms]  Ip[kA] Icd[kA] Ibs[kA] Ioh[kA]      qe    qmin      q0   Vl[V] Te0[keV] Ti0[keV] ne0[e19]   f_ss
          1   0.047    4 5.7e-10       0       1      80       0    1.41    78.5    15.3    5.55    5.55 3.4e+00    0.21    0.21    1.00 1.1e+01
         11    0.25    2 1.8e-12      10       1      80       0     2.4    77.6    15.3    2.89    2.89 2.9e+00    0.30    0.30    1.00 4.1e+00
         21    0.44    2 6.8e-14      20       1      80       0    2.57    77.4    15.3    2.21    2.21 2.7e+00    0.35    0.35    1.00 2.2e+00
         31    0.63    2 1.5e-14      30       1      80       0    2.67    77.3    15.3    1.93    1.93 2.5e+00    0.38    0.38    1.00 1.3e+00
         41    0.81    2 1.7e-14      40       1      80       0    2.72    77.2    15.3    1.79    1.79 2.4e+00    0.39    0.39    1.00 7.5e-01
         51       1    2 2.2e-14      50       1      80       0    2.75    77.2    15.3    1.72    1.72 2.4e+00    0.41    0.41    1.00 4.5e-01

Finally we plot some output quantities

    subplot(121); plot(out.time,out.te0); xlabel('t[s]');ylabel('T_{e,0} [eV]');
    subplot(122); plot(out.rho,out.q); xlabel('\rho_{tor,N}');ylabel('q [-]');