--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/Discovery/Src/vpm.c	Sat Apr 28 11:52:34 2018 +0200
@@ -0,0 +1,2523 @@
+///////////////////////////////////////////////////////////////////////////////
+/// -*- coding: UTF-8 -*-
+///
+/// \file   Discovery/Src/vpm.c
+/// \brief  critical_volume comment by hw
+/// \author Heinrichs Weikamp, Erik C. Baker
+/// \date   19-April-2014
+///
+/// \details
+///
+/// $Id$
+///////////////////////////////////////////////////////////////////////////////
+/// \par Copyright (c) 2014-2018 Heinrichs Weikamp gmbh
+///
+///     This program is free software: you can redistribute it and/or modify
+///     it under the terms of the GNU General Public License as published by
+///     the Free Software Foundation, either version 3 of the License, or
+///     (at your option) any later version.
+///
+///     This program is distributed in the hope that it will be useful,
+///     but WITHOUT ANY WARRANTY; without even the implied warranty of
+///     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+///     GNU General Public License for more details.
+///
+///     You should have received a copy of the GNU General Public License
+///     along with this program.  If not, see <http://www.gnu.org/licenses/>.
+//////////////////////////////////////////////////////////////////////////////
+/// \par Varying Permeability Model (VPM) Decompression Program in c (converted from FORTRAN)
+///
+///     Author:  Erik C. Baker
+///
+///     "DISTRIBUTE FREELY - CREDIT THE AUTHORS"
+///
+///     This program extends the 1986 VPM algorithm (Yount & Hoffman) to include
+///     mixed gas, repetitive, and altitude diving.  Developments to the algorithm
+///     were made by David E. Yount, Eric B. Maiken, and Erik C. Baker over a
+///     period from 1999 to 2001.  This work is dedicated in remembrance of
+///     Professor David E. Yount who passed away on April 27, 2000.
+///
+///     Notes:
+///     1.  This program uses the sixteen (16) half-time compartments of the
+///         Buhlmann ZH-L16 model.  The optional Compartment 1b is used here with
+///         half-times of 1.88 minutes for helium and 5.0 minutes for nitrogen.
+///
+///     2.  This program uses various DEC, IBM, and Microsoft extensions which
+///         may not be supported by all FORTRAN compilers.  Comments are made with
+///         a capital "C" in the first column or an exclamation point "!" placed
+///         in a line after code.  An asterisk "*" in column 6 is a continuation
+///         of the previous line.  All code, except for line numbers, starts in
+///         column 7.
+///
+///     3.  Comments and suggestions for improvements are welcome.  Please
+///         respond by e-mail to:  EBaker@se.aeieng.com
+///
+///     Acknowledgment:  Thanks to Kurt Spaugh for recommendations on how to clean
+///     up the code.
+/// ===============================================================================
+///     Converted to vpmdeco.c using f2c; R.McGinnis (CABER Swe) 5/01
+/// ===============================================================================
+///
+/// ************************ Heirichs Weipkamp **************************************
+///
+/// The original Yount & Baker code has been adjusted for real life calculation.
+///
+/// 1) The original main function has been split in several functions
+///
+/// 2) When the deco zone is reached (while ascending) the gradient factors are kept fix
+///     and critical volume algorithm is switched of. maxfirststopdepth is kept fix
+///     to make shure Boeyls Law algorithm works correctly
+///
+/// 4) gas_loadings_ascent_descend heeds all gaschanges and CCR support has been added
+///
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <math.h>
+#include <time.h>
+
+//#include "compiler.h"
+//#include "sdramc.h"
+#include "vpm.h"
+//#include "buehlmann.h"
+
+#include "decom.h"
+
+//#include "decompression.h"
+//#include "taskmanagement\tissue_calls.h"
+#define true 1
+#define false 0
+
+#define GAS_N2 0
+#define GAS_HE 1
+/* temp vars to simplify UNFMTLISTs */
+float    fO2, fHe, fN2;
+float    dc, rc, ssc;
+short mc;
+
+const _Bool buehlmannSafety = true;
+/* Common Block Declarations */
+
+extern const float SURFACE_TENSION_GAMMA;				//!Adj. Range: 0.015 to 0.065 N/m
+extern const float SKIN_COMPRESSION_GAMMAC;				//!Adj. Range: 0.160 to 0.290 N/m
+extern const float UNITS_FACTOR;
+extern const float WATER_VAPOR_PRESSURE;				// (Schreiner value)  based on respiratory quotien
+extern const float CRIT_VOLUME_PARAMETER_LAMBDA;			//!Adj. Range: 6500 to 8300 fsw-min
+extern const float GRADIENT_ONSET_OF_IMPERM_ATM;			//!Adj. Range: 5.0 to 10.0 atm
+extern const float REGENERATION_TIME_CONSTANT;			//!Adj. Range: 10080 to 51840 min
+extern const float PRESSURE_OTHER_GASES_MMHG;				//!Constant value for PO2 up to 2 atm
+extern const float CONSTANT_PRESSURE_OTHER_GASES; // PRESSURE_OTHER_GASES_MMHG / 760. * UNITS_FACTOR;
+
+extern const float HELIUM_TIME_CONSTANT[];
+extern const float NITROGEN_TIME_CONSTANT[];
+
+float    minimum_deco_stop_time;
+float    run_time, run_time_first_stop;
+float    segment_time;
+short mix_number;
+float    barometric_pressure;
+_Bool altitude_dive_algorithm_off;
+_Bool units_equal_fsw, units_equal_msw;
+
+/* by hw 11.06.2015 to allow */
+float gCNS_VPM;
+
+float    helium_pressure[16], nitrogen_pressure[16];
+//float    helium_pressure_crush[16], nitrogen_pressure_crush[16];
+//float    fraction_helium[MAX_GASMIXES + EXTRA_GASMIXES], fraction_nitrogen[MAX_GASMIXES + EXTRA_GASMIXES];
+//float    initial_critical_radius_he[16], initial_critical_radius_n2[16];
+//float    adjusted_critical_radius_he[16],
+//adjusted_critical_radius_n2[16];
+//float    max_crushing_pressure_he[16],
+//max_crushing_pressure_n2[16];
+float    surface_phase_volume_time[16];
+//float    max_actual_gradient[16];
+//float    amb_pressure_onset_of_imperm[16],
+//gas_tension_onset_of_imperm[16];
+//float    initial_helium_pressure_global[16],
+//initial_nitrogen_pressure_global[16];
+float    regenerated_radius_he[16],
+regenerated_radius_n2[16];
+//float    adjusted_crushing_pressure_he[16],
+//adjusted_crushing_pressure_n2[16];
+float    allowable_gradient_he[16],
+allowable_gradient_n2[16];
+//float    initial_allowable_gradient_he[16],
+//initial_allowable_gradient_n2[16];
+
+//_Bool deco_zone_reached;
+_Bool critical_volume_algorithm_off;
+float max_first_stop_depth;
+float max_deco_ceiling_depth;
+long vpm_time_calc_begin = 0;
+//Boylslaw compensation
+float deco_gradient_he[16];
+float deco_gradient_n2[16];
+int last_nullzeit;
+int vpm_calc_what;
+int count_critical_volume_iteration;
+short number_of_changes;
+ float   depth_change[11];
+ float  step_size_change[11];
+ float    rate_change[11];
+ short mix_change[11];
+
+const _Bool vpm_b = true;
+
+//extern
+/*extern float tissue_N2_saturation[4][16];
+extern float tissue_He_saturation[4][16];
+extern long dv_divetime;
+extern dive_data_t dive_data;
+extern int dive_ambient_pressure_mbar;
+extern long dv_seconds_since_last_dive;
+extern int tts[4];*/
+extern const float float_buehlmann_N2_factor_expositon_20_seconds[];
+extern const float float_buehlmann_He_factor_expositon_20_seconds[];
+extern const float float_buehlmann_N2_factor_expositon_one_minute[];
+extern const float float_buehlmann_He_factor_expositon_one_minute[];
+extern const float float_buehlmann_N2_factor_expositon_five_minutes[];
+extern const float float_buehlmann_He_factor_expositon_five_minutes[];
+extern const float float_buehlmann_N2_factor_expositon_one_hour[];
+extern const float float_buehlmann_He_factor_expositon_one_hour[];
+
+//extern buehlmann_configuration_t buehlmann_config;
+
+extern unsigned char CCR_mode; //0x100 // by tissue_calls.c			// uchar
+
+//extern gaschange2_t gaschange_CCR_backup[2][BUEHLMANN_STRUCT_MAX_GASES];
+
+
+float starting_ambient_pressure_global;
+float ending_ambient_pressure_global;
+float    depth_start_of_deco_calc;
+float    depth_start_of_deco_zone;
+float    first_stop_depth;
+float    run_time_start_of_deco_zone;
+
+float r_nint(float *x);
+float r_int(float *x);
+_Bool repetitive_variables_not_valid = false;
+_Bool nullzeit_unter60;
+//enum VPM_CALC_STATUS{CALC_END,CALC_BEGIN,CALC_NULLZEIT };
+int vpm_calc_status;
+_Bool buehlmann_wait_exceeded = false;
+
+SLifeData* pInput = NULL;
+SVpm* pVpm = NULL;
+SDecoinfo* pDecoInfo = NULL;
+SDiveSettings* pDiveSettings = NULL;
+float r_nint(float *x)
+{
+    return( (*x)>=0 ?
+    floorf(*x + 0.5f) : -floorf(0.5f - *x) );
+}
+
+float r_int(float *x)
+{
+    return( (*x>0) ? floorf(*x) : -floorf(- *x) );
+}
+
+/** private functions
+*/
+int onset_of_impermeability(float *starting_ambient_pressure, float *ending_ambient_pressure, float *rate, short *i);
+int radius_root_finder (float *a, float *b, float *c,float *low_bound,  float *high_bound, float *ending_radius);
+int nuclear_regeneration(float *dive_time);// clock_();
+int calc_deco_ceiling(float *deco_ceiling_depth,_Bool fallowablw);
+
+int calc_barometric_pressure(float *altitude);
+//extern /* Subroutine */ int vpm_repetitive_algorithm();
+//extern /* Subroutine */ int gas_loadings_surface_interval();
+int critical_volume(float *deco_phase_volume_time);	    ;
+int calc_start_of_deco_zone(float *starting_depth, float *rate, float *depth_start_of_deco_zone);
+
+int calc_initial_allowable_gradient(void);
+void decompression_stop(float *deco_stop_depth, float *step_size, _Bool final_deco_calculation);
+int gas_loadings_ascent_descen(float* helium_pressure, float* nitrogen_pressure, float starting_depth,float ending_depth, float rate,_Bool check_gas_change);
+
+int calc_surface_phase_volume_time(void);
+int calc_max_actual_gradient(float *deco_stop_depth);
+int projected_ascent(float *starting_depth, float *rate, float *deco_stop_depth, float *step_size);
+void vpm_calc_deco(void);
+int vpm_calc_critcal_volume(_Bool begin,_Bool calc_nulltime);
+int vpm_check_converged(_Bool calc_nulltime);
+int vpm_calc_final_deco(_Bool begin);
+void BOYLES_LAW_COMPENSATION (float* First_Stop_Depth,float * Deco_Stop_Depth,float* Step_Size);
+int vpm_calc_nullzeit(void);
+int vpm_repetitive_algorithm(float *surface_interval_time);
+void  vpm_init_1(void);
+
+void vpm_calc_deco_ceiling(void);
+
+//extern /* Subroutine */ int gas_loadings_constant_depth();
+//extern /* Subroutine */ int vpm_altitude_dive_algorithm();
+//extern /* Subroutine */ int calc_max_actual_gradient(),
+//projected_ascent();
+
+void ______X_X_X___________________________________________________________(void);
+
+//#define ARGGG
+
+void vpm_reset_variables(void)
+{
+    repetitive_variables_not_valid = true;
+}
+
+void vpm_init_1(void)
+{
+     units_equal_msw = true;
+     units_equal_fsw = false;
+     altitude_dive_algorithm_off= true;			//!Options: ON or OFF
+     minimum_deco_stop_time=1.0;					//!Options: float positive number
+     critical_volume_algorithm_off= false;		//!Options: ON or OFF
+     run_time = 0.;
+    //barometric_pressure = dive_data.surface * 10;
+
+    //mix_number = dive_data.selected_gas + 1;
+
+    max_first_stop_depth = 0;
+    max_deco_ceiling_depth = 0;
+    //deco_zone_reached = false;
+    depth_start_of_deco_calc = 0;
+    depth_start_of_deco_zone = 0;
+    first_stop_depth = 0;
+    run_time_start_of_deco_zone = 0;
+
+    gCNS_VPM = 0;
+}
+
+float vpm_get_CNS(void)
+{
+    return gCNS_VPM;
+}
+
+int  vpm_calc(SLifeData* pINPUT,
+              SDiveSettings* pSettings,
+              SVpm* pVPM,
+              SDecoinfo*
+              pDECOINFO,
+              int calc_what)
+{
+    vpm_init_1();
+    //decom_CreateGasChangeList(pSettings, pINPUT);
+    vpm_calc_what = calc_what;
+    /**clear decoInfo*/
+    pDECOINFO->output_time_to_surface_seconds = 0;
+    pDECOINFO->output_ndl_seconds = 0;
+    pDECOINFO->output_ceiling_meter = 0;
+    pDECOINFO->output_relative_gradient = 0;
+    uint8_t tmp_calc_status;
+    for(int i=0;i<DECOINFO_STRUCT_MAX_STOPS;i++)
+    {
+        pDECOINFO->output_stop_length_seconds[i] = 0;
+    }
+
+    if(pINPUT->dive_time_seconds < 10)
+    {
+        vpm_calc_status = CALC_NULLZEIT;
+        return vpm_calc_status;
+    }
+    pVpm = pVPM;
+    pInput = pINPUT;
+    pDecoInfo = pDECOINFO;
+    pDiveSettings = pSettings;
+
+    if(vpm_calc_status == CALC_NULLZEIT)
+    {
+        tmp_calc_status = vpm_calc_nullzeit();
+    }
+    else
+    {
+        tmp_calc_status = CALC_BEGIN;
+    }
+    //Normal Deco calculation
+    if(tmp_calc_status != CALC_NULLZEIT)
+    {
+        max_first_stop_depth =  pVpm->max_first_stop_depth_save;
+        run_time_start_of_deco_zone = pVpm->run_time_start_of_deco_zone_save;
+        depth_start_of_deco_zone = pVpm->depth_start_of_deco_zone_save;
+        for (int i = 0; i < 16; ++i) {
+            helium_pressure[i] = pInput->tissue_helium_bar[i] * 10;
+            nitrogen_pressure[i] = pInput->tissue_nitrogen_bar[i] * 10;
+        }
+        vpm_calc_deco();
+        tmp_calc_status = vpm_calc_critcal_volume(true,false);
+        if(vpm_calc_what == DECOSTOPS)
+        {
+            pVpm->max_first_stop_depth_save = max_first_stop_depth;
+            pVpm->run_time_start_of_deco_zone_save = run_time_start_of_deco_zone;
+            pVpm->depth_start_of_deco_zone_save = depth_start_of_deco_zone;
+        }
+    }
+
+    //Only Decostops not futute stops
+    if(vpm_calc_what == DECOSTOPS)
+        vpm_calc_status = tmp_calc_status;
+    return vpm_calc_status;
+}
+
+void  vpm_saturation_after_ascent(SLifeData* input)
+{
+    int i = 0;
+    for (i = 0; i < 16; ++i) {
+        pInput->tissue_helium_bar[i] = helium_pressure[i] / 10;
+        pInput->tissue_nitrogen_bar[i] = nitrogen_pressure[i] / 10;
+    }
+    pInput->pressure_ambient_bar = pInput->pressure_surface_bar;
+}
+/* =============================================================================== */
+/*     NOTE ABOUT PRESSURE UNITS USED IN CALCULATIONS: */
+/*     It is the convention in decompression calculations to compute all gas */
+/*     loadings, absolute pressures, partial pressures, etc., in the units of */
+/*     depth pressure that you are diving - either feet of seawater (fsw) or */
+/*     meters of seawater (msw).  This program follows that convention with the */
+/*     the exception that all VPM calculations are performed in SI units (by */
+/*     necessity).  Accordingly, there are several conversions back and forth */
+/*     between the diving pressure units and the SI units. */
+/* =============================================================================== */
+/* =============================================================================== */
+/*     FUNCTION SUBPROGRAM FOR GAS LOADING CALCULATIONS - ASCENT AND DESCENT */
+/* =============================================================================== */
+
+
+
+/* =============================================================================== */
+/*     SUBROUTINE GAS_LOADINGS_ASCENT_DESCENT */
+/*     Purpose: This subprogram applies the Schreiner equation to update the */
+/*     gas loadings (partial pressures of helium and nitrogen) in the half-time */
+/*     compartments due to a linear ascent or descent segment at a constant rate. */
+/* =============================================================================== */
+
+int gas_loadings_ascent_descen(float* helium_pressure,
+                               float* nitrogen_pressure,
+                               float starting_depth,
+                               float ending_depth,
+                               float rate,_Bool check_gas_change)
+{
+    short i;
+    float initial_inspired_n2_pressure,
+    initial_inspired_he_pressure, nitrogen_rate,
+    last_run_time,
+    starting_ambient_pressure,
+    ending_ambient_pressure;
+    float initial_helium_pressure[16];
+    float initial_nitrogen_pressure[16];
+    float helium_rate;
+    float fraction_helium_begin;
+    float fraction_helium_end;
+    float fraction_nitrogen_begin;
+    float fraction_nitrogen_end;
+    float ending_depth_tmp = ending_depth;
+    float segment_time_tmp = 0;
+    /* loop */
+    /* =============================================================================== */
+    /*     CALCULATIONS */
+    /* =============================================================================== */
+    segment_time = (ending_depth_tmp - starting_depth) / rate;
+    last_run_time = run_time;
+    run_time = last_run_time + segment_time;
+    do {
+        ending_depth_tmp = ending_depth;
+        if (starting_depth > ending_depth && check_gas_change && number_of_changes > 1)
+        {
+            for (i = 1; i < number_of_changes; ++i)
+            {
+                if (depth_change[i] < starting_depth && depth_change[i] > ending_depth)
+                {
+                    ending_depth_tmp = depth_change[i];
+                    break;
+                }
+            }
+            for (i = 1; i < number_of_changes; ++i)
+            {
+                if (depth_change[i] >= starting_depth)
+                {
+                    mix_number = mix_change[i];
+                }
+            }
+        }
+        segment_time_tmp = (ending_depth_tmp - starting_depth) / rate;
+        ending_ambient_pressure = ending_depth_tmp + barometric_pressure;
+        starting_ambient_pressure = starting_depth + barometric_pressure;
+        decom_get_inert_gases( starting_ambient_pressure / 10,  (&pDiveSettings->decogaslist[mix_number]), &fraction_nitrogen_begin, &fraction_helium_begin );
+        decom_get_inert_gases( ending_ambient_pressure   / 10,  (&pDiveSettings->decogaslist[mix_number]), &fraction_nitrogen_end, &fraction_helium_end );
+
+        initial_inspired_he_pressure =	(starting_ambient_pressure - WATER_VAPOR_PRESSURE) * fraction_helium_begin;
+        initial_inspired_n2_pressure =	(starting_ambient_pressure - WATER_VAPOR_PRESSURE) * fraction_nitrogen_begin;
+        //helium_rate = *rate * fraction_helium[mix_number - 1];
+        helium_rate = ((ending_ambient_pressure  - WATER_VAPOR_PRESSURE)* fraction_helium_end - initial_inspired_he_pressure)/segment_time_tmp;
+        //nitrogen_rate2 = *rate * fraction_nitrogen[mix_number - 1];
+        nitrogen_rate = ((ending_ambient_pressure  - WATER_VAPOR_PRESSURE)* fraction_nitrogen_end - initial_inspired_n2_pressure)/segment_time_tmp;
+
+
+        decom_oxygen_calculate_cns_stage_SchreinerStyle(segment_time_tmp,&pDiveSettings->decogaslist[mix_number],starting_ambient_pressure/10,ending_ambient_pressure/10,&gCNS_VPM);
+        //if(fabs(nitrogen_rate - nitrogen_rate2) > 0.000001)
+            //return -2;
+        for (i = 1; i <= 16; ++i)
+        {
+            initial_helium_pressure[i - 1] = helium_pressure[i - 1];
+            initial_nitrogen_pressure[i - 1] = nitrogen_pressure[i - 1];
+            helium_pressure[i - 1] =
+            schreiner_equation__2(&initial_inspired_he_pressure,
+            &helium_rate,
+            &segment_time,
+            &HELIUM_TIME_CONSTANT[i - 1],
+            &initial_helium_pressure[i - 1]);
+            nitrogen_pressure[i - 1] =
+            schreiner_equation__2(&initial_inspired_n2_pressure,
+            &nitrogen_rate,
+            &segment_time,
+            &NITROGEN_TIME_CONSTANT[i - 1],
+            &initial_nitrogen_pressure[i - 1]);
+
+            //nextround???
+
+        }
+        starting_depth = ending_depth_tmp;
+    } while(ending_depth_tmp > ending_depth);
+
+    return 0;
+} /* gas_loadings_ascent_descen */
+
+float last_phase_volume_time[16];
+float n2_pressure_start_of_deco_zone[16];
+float he_pressure_start_of_deco_zone[16];
+float phase_volume_time[16];
+float n2_pressure_start_of_ascent[16];
+float  he_pressure_start_of_ascent[16];
+float    run_time_start_of_deco_calc;
+float starting_depth;
+float last_run_time;
+float deco_phase_volume_time;
+
+float run_time_start_of_ascent;
+
+float rate;
+float    step_size;
+_Bool vpm_violates_buehlmann;
+
+void  vpm_calc_deco(void)
+{
+    /* System generated locals */
+
+    //float   deepest_possible_stop_depth;
+//	altitude_of_dive,
+    short i;
+    int j = 0;
+
+    // float    rounding_operation;
+
+    /* =============================================================================== */
+    /*     INPUT PARAMETERS TO BE USED FOR STAGED DECOMPRESSION AND SAVE IN ARRAYS. */
+    /*     ASSIGN INITAL PARAMETERS TO BE USED AT START OF ASCENT */
+    /*     The user has the ability to change mix, ascent rate, and step size in any */
+    /*     combination at any depth during the ascent. */
+    /* =============================================================================== */
+
+    run_time = ((float)pInput->dive_time_seconds )/ 60;
+    count_critical_volume_iteration = 0;
+    number_of_changes = 1;
+
+    barometric_pressure = pInput->pressure_surface_bar * 10;
+    depth_change[0] =(pInput->pressure_ambient_bar - pInput->pressure_surface_bar)* 10;
+    mix_change[0] = 0;
+    rate_change[0 ] = -10;// neu 160215 hw, zuvor: -12;
+    step_size_change[0] = 3;
+    vpm_violates_buehlmann = false;
+
+    for (i = 1; i < BUEHLMANN_STRUCT_MAX_GASES; i++)
+    {
+            depth_change[i] = 0;
+            mix_change[i] = 0;
+    }
+    j = 0;
+
+    for (i = 1; i < BUEHLMANN_STRUCT_MAX_GASES; i++)
+    {
+        if(pDiveSettings->decogaslist[i].change_during_ascent_depth_meter_otherwise_zero  >= depth_change[0] + 1)
+            continue;
+
+        if(pDiveSettings->decogaslist[i].change_during_ascent_depth_meter_otherwise_zero <= 0)
+            break;
+
+        j++;
+        number_of_changes ++;
+        depth_change[j] = pDiveSettings->decogaslist[i].change_during_ascent_depth_meter_otherwise_zero ;
+        mix_change[j] = i;
+        rate_change[j] = -10;// neu 160215 hw, zuvor: -12;
+        step_size_change[j] = 3;
+    }
+
+    starting_depth = depth_change[0] ;
+    mix_number = mix_change[0] ;
+    rate = rate_change[0];
+    step_size = step_size_change[0];
+
+    for (i = 0; i < 16; ++i) {
+        he_pressure_start_of_ascent[i ] = 	helium_pressure[i];
+        n2_pressure_start_of_ascent[i] = 	nitrogen_pressure[i];
+    }
+    run_time_start_of_ascent = run_time;
+    if(starting_depth <= depth_start_of_deco_zone && vpm_calc_what == DECOSTOPS)
+    {
+        pVpm->deco_zone_reached = true;
+        depth_start_of_deco_calc = starting_depth;
+        critical_volume_algorithm_off = true;
+    }
+    else
+    {
+        //if(deco_zone_reached)
+        //{
+            pVpm->deco_zone_reached = false;
+            critical_volume_algorithm_off = false;
+            //max_first_stop_depth = 0;
+            //max_first_stop_depth_save = 0;
+        //}
+        /* =============================================================================== */
+        /*     BEGIN PROCESS OF ASCENT AND DECOMPRESSION */
+        /*     First, calculate the regeneration of critical radii that takes place over */
+        /*     the dive time.  The regeneration time constant has a time scale of weeks */
+        /*     so this will have very little impact on dives of normal length, but will */
+        /*     have major impact for saturation dives. */
+        /* =============================================================================== */
+
+        nuclear_regeneration(&run_time);
+
+        /* =============================================================================== */
+        /*     CALCULATE INITIAL ALLOWABLE GRADIENTS FOR ASCENT */
+        /*     This is based on the maximum effective crushing pressure on critical radii */
+        /*     in each compartment achieved during the dive profile. */
+        /* =============================================================================== */
+
+        calc_initial_allowable_gradient();
+
+        /* =============================================================================== */
+        /*     SAVE VARIABLES AT START OF ASCENT (END OF BOTTOM TIME) SINCE THESE WILL */
+        /*     BE USED LATER TO COMPUTE THE FINAL ASCENT PROFILE THAT IS WRITTEN TO THE */
+        /*     OUTPUT FILE. */
+        /*     The VPM uses an iterative process to compute decompression schedules so */
+        /*     there will be more than one pass through the decompression loop. */
+        /* =============================================================================== */
+
+        /* =============================================================================== */
+        /*     CALCULATE THE DEPTH WHERE THE DECOMPRESSION ZONE BEGINS FOR THIS PROFILE */
+        /*     BASED ON THE INITIAL ASCENT PARAMETERS AND WRITE THE DEEPEST POSSIBLE */
+        /*     DECOMPRESSION STOP DEPTH TO THE OUTPUT FILE */
+        /*     Knowing where the decompression zone starts is very important.  Below */
+        /*     that depth there is no possibility for bubble formation because there */
+        /*     will be no supersaturation gradients.  Deco stops should never start */
+        /*     below the deco zone.  The deepest possible stop deco stop depth is */
+        /*     defined as the next "standard" stop depth above the point where the */
+        /*     leading compartment enters the deco zone.  Thus, the program will not */
+        /*     base this calculation on step sizes larger than 10 fsw or 3 msw.  The */
+        /*     deepest possible stop depth is not used in the program, per se, rather */
+        /*     it is information to tell the diver where to start putting on the brakes */
+        /*     during ascent.  This should be prominently displayed by any deco program. */
+        /* =============================================================================== */
+
+        calc_start_of_deco_zone(&starting_depth, &rate, &depth_start_of_deco_zone);
+        /* =============================================================================== */
+        /*     TEMPORARILY ASCEND PROFILE TO THE START OF THE DECOMPRESSION ZONE, SAVE */
+        /*     VARIABLES AT THIS POINT, AND INITIALIZE VARIABLES FOR CRITICAL VOLUME LOOP */
+        /*     The iterative process of the VPM Critical Volume Algorithm will operate */
+        /*     only in the decompression zone since it deals with excess gas volume */
+        /*     released as a result of supersaturation gradients (not possible below the */
+        /*     decompression zone). */
+        /* =============================================================================== */
+        gas_loadings_ascent_descen(helium_pressure,nitrogen_pressure, starting_depth, depth_start_of_deco_zone, rate, true);
+
+        run_time_start_of_deco_zone = run_time;
+        depth_start_of_deco_calc = depth_start_of_deco_zone;
+
+        for (i = 0; i < 16; ++i)
+        {
+            pVpm->max_actual_gradient[i] = 0.;
+        }
+    }
+
+    for (i = 0; i < 16; ++i)
+    {
+        surface_phase_volume_time[i] = 0.;
+        last_phase_volume_time[i] = 0.;
+        he_pressure_start_of_deco_zone[i] =	helium_pressure[i];
+        n2_pressure_start_of_deco_zone[i] =	nitrogen_pressure[i];
+        //pVpm->max_actual_gradient[i] = 0.;
+    }
+    run_time_start_of_deco_calc = run_time;
+}
+/* =============================================================================== */
+/*     START OF CRITICAL VOLUME LOOP */
+/*     This loop operates between Lines 50 and 100.  If the Critical Volume */
+/*     Algorithm is toggled "off" in the program settings, there will only be */
+/*     one pass through this loop.  Otherwise, there will be two or more passes */
+/*     through this loop until the deco schedule is "converged" - that is when a */
+/*     comparison between the phase volume time of the present iteration and the */
+/*     last iteration is less than or equal to one minute.  This implies that */
+/*     the volume of released gas in the most recent iteration differs from the */
+/*     "critical" volume limit by an acceptably small amount.  The critical */
+/*     volume limit is set by the Critical Volume Parameter Lambda in the program */
+/*     settings (default setting is 7500 fsw-min with adjustability range from */
+/*     from 6500 to 8300 fsw-min according to Bruce Wienke). */
+/* =============================================================================== */
+/* L50: */
+
+float  deco_stop_depth;
+int vpm_calc_critcal_volume(_Bool begin,
+                            _Bool calc_nulltime)
+{   /* loop will run continuous there is an exit stateme */
+
+    short i;
+
+    float rounding_operation2;
+    //float    ending_depth;
+    float deco_ceiling_depth;
+
+    //float deco_time;
+    int count = 0;
+    _Bool first_stop;
+    int dp = 0;
+    float tissue_He_saturation[16];
+    float tissue_N2_saturation[16];
+    float vpm_buehlmann_safety_gradient = 1.0f - (((float)pDiveSettings->vpm_conservatism) / 40);
+    /* =============================================================================== */
+    /*     CALCULATE CURRENT DECO CEILING BASED ON ALLOWABLE SUPERSATURATION */
+    /*     GRADIENTS AND SET FIRST DECO STOP.  CHECK TO MAKE SURE THAT SELECTED STEP */
+    /*     SIZE WILL NOT ROUND UP FIRST STOP TO A DEPTH THAT IS BELOW THE DECO ZONE. */
+    /* =============================================================================== */
+    if(begin)
+    {
+        if(depth_start_of_deco_calc < max_first_stop_depth )
+        {
+            if(vpm_b)
+            {
+                BOYLES_LAW_COMPENSATION(&max_first_stop_depth, &depth_start_of_deco_calc, &step_size);
+            }
+            calc_deco_ceiling(&deco_ceiling_depth, false);
+        }
+        else
+            calc_deco_ceiling(&deco_ceiling_depth, true);
+
+
+        if (deco_ceiling_depth <= 0.0f) {
+            deco_stop_depth = 0.0f;
+        } else {
+            rounding_operation2 = deco_ceiling_depth / step_size + ( float)0.5f;
+            deco_stop_depth = r_nint(&rounding_operation2) * step_size;
+        }
+
+        // buehlmann safety
+        if(buehlmannSafety)
+        {
+            for (i = 0; i < 16; i++)
+            {
+                tissue_He_saturation[i] = helium_pressure[i] / 10;
+                tissue_N2_saturation[i] = nitrogen_pressure[i] / 10;
+            }
+
+            if(!decom_tissue_test_tolerance(tissue_N2_saturation, tissue_He_saturation, vpm_buehlmann_safety_gradient, (deco_stop_depth / 10.0f) + pInput->pressure_surface_bar))
+            {
+
+                vpm_violates_buehlmann = true;
+                do {
+                    deco_stop_depth += 3;
+                } while (!decom_tissue_test_tolerance(tissue_N2_saturation, tissue_He_saturation, vpm_buehlmann_safety_gradient, (deco_stop_depth / 10.0f) + pInput->pressure_surface_bar));
+            }
+        }
+
+        /* =============================================================================== */
+        /*     PERFORM A SEPARATE "PROJECTED ASCENT" OUTSIDE OF THE MAIN PROGRAM TO MAKE */
+        /*     SURE THAT AN INCREASE IN GAS LOADINGS DURING ASCENT TO THE FIRST STOP WILL */
+        /*     NOT CAUSE A VIOLATION OF THE DECO CEILING.  IF SO, ADJUST THE FIRST STOP */
+        /*     DEEPER BASED ON STEP SIZE UNTIL A SAFE ASCENT CAN BE MADE. */
+        /*     Note: this situation is a possibility when ascending from extremely deep */
+        /*     dives or due to an unusual gas mix selection. */
+        /*     CHECK AGAIN TO MAKE SURE THAT ADJUSTED FIRST STOP WILL NOT BE BELOW THE */
+        /*     DECO ZONE. */
+        /* =============================================================================== */
+        if (deco_stop_depth < depth_start_of_deco_calc)
+        {
+            projected_ascent(&depth_start_of_deco_calc, &rate, &deco_stop_depth, &step_size);
+        }
+
+        /*if (deco_stop_depth > depth_start_of_deco_zone) {
+            printf("\t\n");
+            printf(fmt_905);
+            printf(fmt_900);
+            printf("\nPROGRAM TERMINATED\n");
+            exit(1);
+        }*/
+
+        /* =============================================================================== */
+        /*     HANDLE THE SPECIAL CASE WHEN NO DECO STOPS ARE REQUIRED - ASCENT CAN BE */
+        /*     MADE DIRECTLY TO THE SURFACE */
+        /*     Write ascent data to output file and exit the Critical Volume Loop. */
+        /* =============================================================================== */
+
+        if (deco_stop_depth == 0.0f)
+        {
+            if(calc_nulltime)
+            {
+                return CALC_END;
+            }
+            if(pVpm->deco_zone_reached)
+            {
+                for(dp = 0;dp < DECOINFO_STRUCT_MAX_STOPS;dp++)
+                {
+                    pDecoInfo->output_stop_length_seconds[dp] = 0;
+                }
+                pDecoInfo->output_ndl_seconds = 0;
+                /*max_first_stop_depth = 0;
+                deco_zone_reached = false;
+                depth_start_of_deco_calc = 0;
+                depth_start_of_deco_zone = 0;
+                first_stop_depth = 0;
+                max_first_stop_depth_save = 0;
+                depth_start_of_deco_zone_save = 0;
+                run_time_start_of_deco_zone_save = 0;
+                tts[DECOSTOPS] = 0;
+                tts[NULLZEIT] = 0;
+                tts[FUTURESTOPS] = 0;
+                nullzeit_unter60 = false;
+                vpm_calc_status = CALC_NULLZEIT;
+                vpm_calc_what = DECOSTOPS;
+                float surfacetime = 0;
+                vpm_repetitive_algorithm(&surfacetime);*/
+
+            }
+
+            return CALC_NULLZEIT;
+                    /* exit the critical volume l */
+        }
+
+    /* =============================================================================== */
+    /*     ASSIGN VARIABLES FOR ASCENT FROM START OF DECO ZONE TO FIRST STOP.  SAVE */
+    /*     FIRST STOP DEPTH FOR LATER USE WHEN COMPUTING THE FINAL ASCENT PROFILE */
+    /* =============================================================================== */
+        deco_stop_depth = fmaxf(deco_stop_depth,(float)pDiveSettings->last_stop_depth_bar * 10);
+        starting_depth = depth_start_of_deco_calc;
+        first_stop_depth = deco_stop_depth;
+        first_stop = true;
+    }
+    /* =============================================================================== */
+    /*     DECO STOP LOOP BLOCK WITHIN CRITICAL VOLUME LOOP */
+    /*     This loop computes a decompression schedule to the surface during each */
+    /*     iteration of the critical volume loop.  No output is written from this */
+    /*     loop, rather it computes a schedule from which the in-water portion of the */
+    /*     total phase volume time (Deco_Phase_Volume_Time) can be extracted.  Also, */
+    /*     the gas loadings computed at the end of this loop are used the subroutine */
+    /*     which computes the out-of-water portion of the total phase volume time */
+    /*     (Surface_Phase_Volume_Time) for that schedule. */
+
+    /*     Note that exit is made from the loop after last ascent is made to a deco */
+    /*     stop depth that is less than or equal to zero.  A final deco stop less */
+    /*     than zero can happen when the user makes an odd step size change during */
+    /*     ascent - such as specifying a 5 msw step size change at the 3 msw stop! */
+    /* =============================================================================== */
+
+    while(true)           /* loop will run continuous there is an break statement */
+    {
+        if(starting_depth > deco_stop_depth )
+            gas_loadings_ascent_descen(helium_pressure, nitrogen_pressure, starting_depth, deco_stop_depth, rate,first_stop);
+
+        first_stop = false;
+        if (deco_stop_depth <= 0.0f)
+        {
+            break;
+        }
+        if (number_of_changes > 1)
+        {
+            int i1 = number_of_changes;
+            for (i = 2; i <= i1; ++i) {
+                if (depth_change[i - 1] >= deco_stop_depth)
+                {
+                    mix_number = mix_change[i - 1];
+                    rate = rate_change[i - 1];
+                    step_size = step_size_change[i - 1];
+                }
+            }
+        }
+        if(vpm_b)
+        {
+            float fist_stop_depth2 =  fmaxf(first_stop_depth,max_first_stop_depth);
+            BOYLES_LAW_COMPENSATION(&fist_stop_depth2, &deco_stop_depth, &step_size);
+        }
+        decompression_stop(&deco_stop_depth, &step_size, false);
+        starting_depth = deco_stop_depth;
+
+        if(deco_stop_depth == (float)pDiveSettings->last_stop_depth_bar * 10)
+            deco_stop_depth = 0;
+        else
+        {
+            deco_stop_depth = deco_stop_depth - step_size;
+            deco_stop_depth = fmaxf(deco_stop_depth,(float)pDiveSettings->last_stop_depth_bar * 10);
+        }
+
+        count++;
+        //if(count > 14)
+            //return CALC_CRITICAL2;
+        /* L60: */
+    }
+
+        return vpm_check_converged(calc_nulltime);
+}
+/* =============================================================================== */
+/*     COMPUTE TOTAL PHASE VOLUME TIME AND MAKE CRITICAL VOLUME COMPARISON */
+/*     The deco phase volume time is computed from the run time.  The surface */
+/*     phase volume time is computed in a subroutine based on the surfacing gas */
+/*     loadings from previous deco loop block.  Next the total phase volume time */
+/*     (in-water + surface) for each compartment is compared against the previous */
+/*     total phase volume time.  The schedule is converged when the difference is */
+/*     less than or equal to 1 minute in any one of the 16 compartments. */
+
+/*     Note:  the "phase volume time" is somewhat of a mathematical concept. */
+/*     It is the time divided out of a total integration of supersaturation */
+/*     gradient x time (in-water and surface).  This integration is multiplied */
+/*     by the excess bubble number to represent the amount of free-gas released */
+/*     as a result of allowing a certain number of excess bubbles to form. */
+/* =============================================================================== */
+/* end of deco stop loop */
+
+int vpm_check_converged(_Bool calc_nulltime)
+{
+
+    short i;
+    float    critical_volume_comparison;
+    float    r1;
+    _Bool schedule_converged = false;
+
+
+    deco_phase_volume_time = run_time - run_time_start_of_deco_zone;
+    calc_surface_phase_volume_time();
+
+    for (i = 1; i <= 16; ++i)
+    {
+        phase_volume_time[i - 1] =
+        deco_phase_volume_time + surface_phase_volume_time[i - 1];
+        critical_volume_comparison = (r1 = phase_volume_time[i - 1] - last_phase_volume_time[i - 1], fabs(r1));
+
+        if (critical_volume_comparison <= 1.0f)
+        {
+            schedule_converged = true;
+        }
+    }
+
+/* =============================================================================== */
+/*     CRITICAL VOLUME DECISION TREE BETWEEN LINES 70 AND 99 */
+/*     There are two options here.  If the Critical Volume Agorithm setting is */
+/*     "on" and the schedule is converged, or the Critical Volume Algorithm */
+/*     setting was "off" in the first place, the program will re-assign variables */
+/*     to their values at the start of ascent (end of bottom time) and process */
+/*     a complete decompression schedule once again using all the same ascent */
+/*     parameters and first stop depth.  This decompression schedule will match */
+/*     the last iteration of the Critical Volume Loop and the program will write */
+/*     the final deco schedule to the output file. */
+
+/*     Note: if the Critical Volume Agorithm setting was "off", the final deco */
+/*     schedule will be based on "Initial Allowable Supersaturation Gradients." */
+/*     If it was "on", the final schedule will be based on "Adjusted Allowable */
+/*     Supersaturation Gradients" (gradients that are "relaxed" as a result of */
+/*     the Critical Volume Algorithm). */
+
+/*     If the Critical Volume Agorithm setting is "on" and the schedule is not */
+/*     converged, the program will re-assign variables to their values at the */
+/*     start of the deco zone and process another trial decompression schedule. */
+/* =============================================================================== */
+/* L70: */
+    //Not more than 4 iteration allowed
+    count_critical_volume_iteration++;
+    if(count_critical_volume_iteration > 4)
+    {
+        //return CALC_FINAL_DECO;
+        if(calc_nulltime)
+            return CALC_FINAL_DECO;
+        else
+            return vpm_calc_final_deco(true);
+    }
+    if (schedule_converged || critical_volume_algorithm_off)
+    {
+
+        //return CALC_FINAL_DECO;
+        if(calc_nulltime)
+            return CALC_FINAL_DECO;
+        else
+            return vpm_calc_final_deco(true);
+/* final deco schedule */
+/* exit critical volume l */
+
+/* =============================================================================== */
+/*     IF SCHEDULE NOT CONVERGED, COMPUTE RELAXED ALLOWABLE SUPERSATURATION */
+/*     GRADIENTS WITH VPM CRITICAL VOLUME ALGORITHM AND PROCESS ANOTHER */
+/*     ITERATION OF THE CRITICAL VOLUME LOOP */
+/* =============================================================================== */
+
+    } else {
+        critical_volume(&deco_phase_volume_time);
+        deco_phase_volume_time = 0.;
+        run_time = run_time_start_of_deco_calc;
+        starting_depth = depth_start_of_deco_calc;
+        mix_number = mix_change[0];
+        rate = rate_change[0];
+        step_size = step_size_change[0];
+        for (i = 1; i <= 16; ++i)
+        {
+                last_phase_volume_time[i - 1] = phase_volume_time[i - 1];
+                helium_pressure[i - 1] = he_pressure_start_of_deco_zone[i - 1];
+                nitrogen_pressure[i - 1] = n2_pressure_start_of_deco_zone[i - 1];
+        }
+        if(calc_nulltime)
+            return CALC_CRITICAL;
+        else
+            return vpm_calc_critcal_volume(true, false);
+    }
+/* end of critical volume decision */
+/* L100: */
+    //  }/* end of critical vol loop */
+}
+
+void vpm_calc_deco_ceiling(void)
+{
+
+    short i;
+// hw 1601209	 float    r1;
+// hw 1601209	 float    stop_time;
+// hw 1601209	 int count = 0;
+    //static int dp_max;
+    //static float surfacetime;
+    // _Bool first_stop = false;
+    float tissue_He_saturation[16];
+    float tissue_N2_saturation[16];
+    float vpm_buehlmann_safety_gradient = 1.0f - (((float)pDiveSettings->vpm_conservatism) / 40);
+    //max_first_stop_depth = fmaxf(first_stop_depth,max_first_stop_depth);
+
+    /** CALC DECO Ceiling ******************************************************************/
+    /** Not when Future stops */
+    if(vpm_calc_what == DECOSTOPS)
+    {
+
+        for (i = 1; i <= 16; ++i)
+        {
+            helium_pressure[i - 1] = he_pressure_start_of_deco_zone[i - 1];
+            nitrogen_pressure[i - 1] = n2_pressure_start_of_deco_zone[i - 1];
+        }
+        run_time = run_time_start_of_ascent;// run_time_start_of_ascent;
+        starting_depth = depth_change[0];
+        mix_number = mix_change[0];
+        rate = rate_change[0];
+        //gas_loadings_ascent_descen(helium_pressure,nitrogen_pressure, starting_depth, depth_start_of_deco_calc, rate, true);
+
+        float deco_ceiling_depth = 0.0f;
+        if(depth_start_of_deco_calc > max_deco_ceiling_depth)
+        {
+            calc_deco_ceiling(&deco_ceiling_depth, true);
+        }
+        if(buehlmannSafety)
+        {
+            for (i = 0; i < 16; i++)
+            {
+                tissue_He_saturation[i] = helium_pressure[i] / 10;
+                tissue_N2_saturation[i] = nitrogen_pressure[i] / 10;
+            }
+
+            if(!decom_tissue_test_tolerance(tissue_N2_saturation, tissue_He_saturation, vpm_buehlmann_safety_gradient, (deco_ceiling_depth / 10.0f) + pInput->pressure_surface_bar))
+            {
+
+                vpm_violates_buehlmann = true;
+                do {
+                    deco_ceiling_depth += 0.1f;
+                } while (!decom_tissue_test_tolerance(tissue_N2_saturation, tissue_He_saturation, vpm_buehlmann_safety_gradient, (deco_ceiling_depth / 10.0f) + pInput->pressure_surface_bar));
+            }
+        }
+
+        if (deco_ceiling_depth < depth_start_of_deco_calc)
+        {
+            projected_ascent(&depth_start_of_deco_calc, &rate, &deco_ceiling_depth, &step_size);
+        }
+
+        max_deco_ceiling_depth = fmaxf(max_deco_ceiling_depth,deco_ceiling_depth);
+
+        if(depth_start_of_deco_calc > deco_ceiling_depth)
+        {
+            gas_loadings_ascent_descen(helium_pressure,nitrogen_pressure, depth_start_of_deco_calc,deco_ceiling_depth, rate, true);
+            //surfacetime += segment_time;
+        }
+
+        if(vpm_b)
+        {
+            BOYLES_LAW_COMPENSATION(&max_deco_ceiling_depth, &deco_ceiling_depth, &step_size);
+        }
+        calc_deco_ceiling(&deco_ceiling_depth, false);
+
+        // buehlmann safety
+        if(vpm_violates_buehlmann)
+        {
+            for (i = 0; i < 16; i++)
+            {
+                tissue_He_saturation[i] = helium_pressure[i] / 10;
+                tissue_N2_saturation[i] = nitrogen_pressure[i] / 10;
+            }
+
+            if(!decom_tissue_test_tolerance(tissue_N2_saturation, tissue_He_saturation, vpm_buehlmann_safety_gradient, (deco_ceiling_depth / 10.0f) + pInput->pressure_surface_bar))
+            {
+
+                vpm_violates_buehlmann = true;
+                do {
+                    deco_ceiling_depth += 0.1f;
+                } while (!decom_tissue_test_tolerance(tissue_N2_saturation, tissue_He_saturation, vpm_buehlmann_safety_gradient, (deco_ceiling_depth / 10.0f) + pInput->pressure_surface_bar));
+            }
+        }
+        // output_ceiling_meter
+        if(deco_ceiling_depth > first_stop_depth)
+            deco_ceiling_depth = first_stop_depth;
+        pDecoInfo->output_ceiling_meter = deco_ceiling_depth ;
+    }
+    else
+    {
+         pDecoInfo->output_ceiling_meter = 0;
+    }
+
+    // fix hw 160627
+    if(pDecoInfo->output_ceiling_meter < 0)
+        pDecoInfo->output_ceiling_meter = 0;
+
+    /*** End CALC ceiling    ***************************************************/
+}
+
+
+/* =============================================================================== */
+/*     DECO STOP LOOP BLOCK FOR FINAL DECOMPRESSION SCHEDULE */
+/* =============================================================================== */
+
+int vpm_calc_final_deco(_Bool begin)
+{
+    short i;
+    float    r1;
+    float    stop_time;
+    int count = 0;
+    static int dp_max;
+    static float surfacetime;
+    _Bool first_stop = false;
+// hw 1601209 float tissue_He_saturation[16];
+// hw 1601209	float tissue_N2_saturation[16];
+    float vpm_buehlmann_safety_gradient = 1.0f - (((float)pDiveSettings->vpm_conservatism) / 40);
+    max_first_stop_depth = fmaxf(first_stop_depth,max_first_stop_depth);
+    if(begin)
+    {
+        gCNS_VPM = 0;
+        dp_max = 0;
+        for (i = 1; i <= 16; ++i)
+        {
+            helium_pressure[i - 1] =
+            he_pressure_start_of_ascent[i - 1];
+            nitrogen_pressure[i - 1] =
+            n2_pressure_start_of_ascent[i - 1];
+        }
+        run_time = run_time_start_of_ascent;// run_time_start_of_ascent;
+        starting_depth = depth_change[0];
+        mix_number = mix_change[0];
+        rate = rate_change[0];
+        step_size = step_size_change[0];
+        deco_stop_depth = first_stop_depth;
+        max_first_stop_depth = fmaxf(first_stop_depth,max_first_stop_depth);
+        last_run_time = 0.;
+
+
+
+        /* =============================================================================== */
+        /*     DECO STOP LOOP BLOCK FOR FINAL DECOMPRESSION SCHEDULE */
+        /* =============================================================================== */
+        surfacetime = 0;
+        first_stop = true;
+    }
+
+    while(true)       /* loop will run continuous until there is an break statement */
+    {
+        if(starting_depth > deco_stop_depth)
+        {
+            gas_loadings_ascent_descen(helium_pressure,nitrogen_pressure, starting_depth,deco_stop_depth, rate, first_stop);
+            surfacetime += segment_time;
+        }
+
+        /* =============================================================================== */
+        /*     DURING FINAL DECOMPRESSION SCHEDULE PROCESS, COMPUTE MAXIMUM ACTUAL */
+        /*     SUPERSATURATION GRADIENT RESULTING IN EACH COMPARTMENT */
+        /*     If there is a repetitive dive, this will be used later in the VPM */
+        /*     Repetitive Algorithm to adjust the values for critical radii. */
+        /* =============================================================================== */
+        if(vpm_calc_what == DECOSTOPS)
+            calc_max_actual_gradient(&deco_stop_depth);
+
+        if (deco_stop_depth <= 0.0f) {
+            break;
+        }
+        if (number_of_changes > 1)
+        {
+            int i1 = number_of_changes;
+            for (i = 2; i <= i1; ++i)
+            {
+                if (depth_change[i - 1] >= deco_stop_depth)
+                {
+                    mix_number = mix_change[i - 1];
+                    rate = rate_change[i - 1];
+                    step_size = step_size_change[i - 1];
+                }
+            }
+        }
+
+        if(first_stop)
+        {
+            run_time_first_stop = run_time;
+            first_stop = false;
+        }
+        if(vpm_b)
+        {
+            BOYLES_LAW_COMPENSATION(&max_first_stop_depth, &deco_stop_depth, &step_size);
+        }
+        decompression_stop(&deco_stop_depth, &step_size, true);
+
+        /* =============================================================================== */
+        /*     This next bit justs rounds up the stop time at the first stop to be in */
+        /*     whole increments of the minimum stop time (to make for a nice deco table). */
+        /* =============================================================================== */
+
+        if (last_run_time == 0.0f)
+        {
+            r1 = segment_time / minimum_deco_stop_time + 0.5f;
+            stop_time = r_int(&r1) * minimum_deco_stop_time;
+        } else {
+            stop_time = run_time - last_run_time;
+        }
+        stop_time = segment_time;
+        surfacetime += stop_time;
+        if((vpm_calc_what == DECOSTOPS) || (vpm_calc_what == BAILOUTSTOPS))
+        {
+            int dp = 0;
+            if(deco_stop_depth == (float)pDiveSettings->last_stop_depth_bar * 10)
+            {
+                dp = 0;
+            }
+            else
+            {
+                    dp = 1 + (int)((deco_stop_depth - (pDiveSettings->input_second_to_last_stop_depth_bar * 10)) / step_size);
+            }
+            dp_max = (int)fmaxf(dp_max,dp);
+            if(dp < DECOINFO_STRUCT_MAX_STOPS)
+            {
+                int stop_time_seconds = fminf((999 * 60), (int)(stop_time *60));
+                //
+
+                //if(vpm_calc_what == DECOSTOPS)
+                    pDecoInfo->output_stop_length_seconds[dp] = (unsigned short)stop_time_seconds;
+                //else
+                    //decostop_bailout[dp] = (unsigned short)stop_time_seconds;
+            }
+        }
+
+
+        /* =============================================================================== */
+        /*     DURING FINAL DECOMPRESSION SCHEDULE, IF MINIMUM STOP TIME PARAMETER IS A */
+        /*     WHOLE NUMBER (i.e. 1 minute) THEN WRITE DECO SCHEDULE USING short */
+        /*     NUMBERS (looks nicer).  OTHERWISE, USE DECIMAL NUMBERS. */
+        /*     Note: per the request of a noted exploration diver(!), program now allows */
+        /*     a minimum stop time of less than one minute so that total ascent time can */
+        /*     be minimized on very long dives.  In fact, with step size set at 1 fsw or */
+        /*     0.2 msw and minimum stop time set at 0.1 minute (6 seconds), a near */
+        /*     continuous decompression schedule can be computed. */
+        /* =============================================================================== */
+
+        starting_depth = deco_stop_depth;
+        if(deco_stop_depth == (float)pDiveSettings->last_stop_depth_bar * 10)
+            deco_stop_depth = 0;
+        else
+        {
+            deco_stop_depth = deco_stop_depth - step_size;
+            deco_stop_depth = fmaxf(deco_stop_depth,(float)pDiveSettings->last_stop_depth_bar * 10);
+        }
+
+        last_run_time = run_time;
+        count++;
+        //if(count > 14)
+            //return CALC_FINAL_DECO2;
+        /* L80: */
+    } /* for final deco sche */
+
+    if( (vpm_calc_what == DECOSTOPS)  || (vpm_calc_what == BAILOUTSTOPS))
+    {
+        for(int dp = dp_max +1;dp < DECOINFO_STRUCT_MAX_STOPS;dp++)
+        {
+            //if(vpm_calc_what == DECOSTOPS)
+                pDecoInfo->output_stop_length_seconds[dp] = 0;
+            //else
+                //decostop_bailout[dp] = 0;
+        }
+    }
+    pDecoInfo->output_time_to_surface_seconds = (int)(surfacetime * 60);
+    pDecoInfo->output_ndl_seconds = 0;
+
+    vpm_calc_deco_ceiling();
+    /* end of deco stop lo */
+    return CALC_END;
+}
+
+/* =============================================================================== */
+/*     SUBROUTINE NUCLEAR_REGENERATION */
+/*     Purpose: This subprogram calculates the regeneration of VPM critical */
+/*     radii that takes place over the dive time.  The regeneration time constant */
+/*     has a time scale of weeks so this will have very little impact on dives of */
+/*     normal length, but will have a major impact for saturation dives. */
+/* =============================================================================== */
+
+int nuclear_regeneration(float *dive_time)
+{
+    /* Local variables */
+    float crush_pressure_adjust_ratio_he,
+    ending_radius_n2,
+    ending_radius_he;
+    short i;
+    float crushing_pressure_pascals_n2,
+    crushing_pressure_pascals_he,
+    adj_crush_pressure_n2_pascals,
+    adj_crush_pressure_he_pascals,
+    crush_pressure_adjust_ratio_n2;
+
+    /* loop */
+    /* =============================================================================== */
+    /*     CALCULATIONS */
+    /*     First convert the maximum crushing pressure obtained for each compartment */
+    /*     to Pascals.  Next, compute the ending radius for helium and nitrogen */
+    /*     critical nuclei in each compartment. */
+    /* =============================================================================== */
+
+    for (i = 1; i <= 16; ++i)
+    {
+        crushing_pressure_pascals_he =
+        pVpm->max_crushing_pressure_he[i - 1] / UNITS_FACTOR * 101325.0f;
+        crushing_pressure_pascals_n2 =
+        pVpm->max_crushing_pressure_n2[i - 1] / UNITS_FACTOR * 101325.0f;
+        ending_radius_he =
+        1.0f / (crushing_pressure_pascals_he /
+        ((SKIN_COMPRESSION_GAMMAC - SURFACE_TENSION_GAMMA) * 2.0f) +
+        1.0f / pVpm->adjusted_critical_radius_he[i - 1]);
+        ending_radius_n2 =
+        1.0f / (crushing_pressure_pascals_n2 /
+        ((SKIN_COMPRESSION_GAMMAC - SURFACE_TENSION_GAMMA) * 2.0f) +
+        1.0f / pVpm->adjusted_critical_radius_n2[i - 1]);
+
+        /* =============================================================================== */
+        /*     A "regenerated" radius for each nucleus is now calculated based on the */
+        /*     regeneration time constant.  This means that after application of */
+        /*     crushing pressure and reduction in radius, a nucleus will slowly grow */
+        /*     back to its original initial radius over a period of time.  This */
+        /*     phenomenon is probabilistic in nature and depends on absolute temperature. */
+        /*     It is independent of crushing pressure. */
+        /* =============================================================================== */
+
+        regenerated_radius_he[i - 1] =
+        pVpm->adjusted_critical_radius_he[i - 1] +
+        (ending_radius_he - pVpm->adjusted_critical_radius_he[i - 1]) *
+        expf(-(*dive_time) / REGENERATION_TIME_CONSTANT);
+        regenerated_radius_n2[i - 1] =
+        pVpm->adjusted_critical_radius_n2[i - 1] +
+        (ending_radius_n2 - pVpm->adjusted_critical_radius_n2[i - 1]) *
+        expf(-(*dive_time) / REGENERATION_TIME_CONSTANT);
+
+        /* =============================================================================== */
+        /*     In order to preserve reference back to the initial critical radii after */
+        /*     regeneration, an "adjusted crushing pressure" for the nuclei in each */
+        /*     compartment must be computed.  In other words, this is the value of */
+        /*     crushing pressure that would have reduced the original nucleus to the */
+        /*     to the present radius had regeneration not taken place.  The ratio */
+        /*     for adjusting crushing pressure is obtained from algebraic manipulation */
+        /*     of the standard VPM equations.  The adjusted crushing pressure, in lieu */
+        /*     of the original crushing pressure, is then applied in the VPM Critical */
+        /*     Volume Algorithm and the VPM Repetitive Algorithm. */
+        /* =============================================================================== */
+
+        crush_pressure_adjust_ratio_he =
+        ending_radius_he * (pVpm->adjusted_critical_radius_he[i - 1] -
+        regenerated_radius_he[i - 1]) /
+        (regenerated_radius_he[i - 1] *
+        (pVpm->adjusted_critical_radius_he[i - 1] -
+        ending_radius_he));
+        crush_pressure_adjust_ratio_n2 =
+        ending_radius_n2 * (pVpm->adjusted_critical_radius_n2[i - 1] -
+        regenerated_radius_n2[i - 1]) /
+        (regenerated_radius_n2[i - 1] *
+        (pVpm->adjusted_critical_radius_n2[i - 1] -
+        ending_radius_n2));
+        adj_crush_pressure_he_pascals =
+        crushing_pressure_pascals_he * crush_pressure_adjust_ratio_he;
+        adj_crush_pressure_n2_pascals =
+        crushing_pressure_pascals_n2 * crush_pressure_adjust_ratio_n2;
+        pVpm->adjusted_crushing_pressure_he[i - 1] =
+        adj_crush_pressure_he_pascals / 101325.0f * UNITS_FACTOR;
+        pVpm->adjusted_crushing_pressure_n2[i - 1] =
+        adj_crush_pressure_n2_pascals / 101325.0f * UNITS_FACTOR;
+    }
+    return 0;
+} /* nuclear_regeneration */
+
+/* =============================================================================== */
+/*     SUBROUTINE CALC_INITIAL_ALLOWABLE_GRADIENT */
+/*     Purpose: This subprogram calculates the initial allowable gradients for */
+/*     helium and nitrogren in each compartment.  These are the gradients that */
+/*     will be used to set the deco ceiling on the first pass through the deco */
+/*     loop.  If the Critical Volume Algorithm is set to "off", then these */
+/*     gradients will determine the final deco schedule.  Otherwise, if the */
+/*     Critical Volume Algorithm is set to "on", these gradients will be further */
+/*     "relaxed" by the Critical Volume Algorithm subroutine.  The initial */
+/*     allowable gradients are referred to as "PssMin" in the papers by Yount */
+/*     and colleauges, i.e., the minimum supersaturation pressure gradients */
+/*     that will probe bubble formation in the VPM nuclei that started with the */
+/*     designated minimum initial radius (critical radius). */
+
+/*     The initial allowable gradients are computed directly from the */
+/*     "regenerated" radii after the Nuclear Regeneration subroutine.  These */
+/*     gradients are tracked separately for helium and nitrogen. */
+/* =============================================================================== */
+
+int calc_initial_allowable_gradient()
+{
+    float initial_allowable_grad_n2_pa,
+    initial_allowable_grad_he_pa;
+    short i;
+
+    /* loop */
+    /* =============================================================================== */
+    /*     CALCULATIONS */
+    /*     The initial allowable gradients are computed in Pascals and then converted */
+    /*     to the diving pressure units.  Two different sets of arrays are used to */
+    /*     save the calculations - Initial Allowable Gradients and Allowable */
+    /*     Gradients.  The Allowable Gradients are assigned the values from Initial */
+    /*     Allowable Gradients however the Allowable Gradients can be changed later */
+    /*     by the Critical Volume subroutine.  The values for the Initial Allowable */
+    /*     Gradients are saved in a global array for later use by both the Critical */
+    /*     Volume subroutine and the VPM Repetitive Algorithm subroutine. */
+    /* =============================================================================== */
+
+    for (i = 1; i <= 16; ++i)
+    {
+        initial_allowable_grad_n2_pa =
+            SURFACE_TENSION_GAMMA * 2.0f *
+            (SKIN_COMPRESSION_GAMMAC - SURFACE_TENSION_GAMMA) /
+            (regenerated_radius_n2[i - 1] * SKIN_COMPRESSION_GAMMAC);
+
+        initial_allowable_grad_he_pa =
+            SURFACE_TENSION_GAMMA * 2.0f *
+            (SKIN_COMPRESSION_GAMMAC - SURFACE_TENSION_GAMMA) /
+            (regenerated_radius_he[i - 1] * SKIN_COMPRESSION_GAMMAC);
+
+        pVpm->initial_allowable_gradient_n2[i - 1] =
+            initial_allowable_grad_n2_pa / 101325.0f * UNITS_FACTOR;
+
+        pVpm->initial_allowable_gradient_he[i - 1] =
+            initial_allowable_grad_he_pa / 101325.0f * UNITS_FACTOR;
+
+        allowable_gradient_he[i - 1] =
+            pVpm->initial_allowable_gradient_he[i - 1];
+
+        allowable_gradient_n2[i - 1] =
+            pVpm->initial_allowable_gradient_n2[i - 1];
+    }
+    return 0;
+} /* calc_initial_allowable_gradient */
+
+/* =============================================================================== */
+/*     SUBROUTINE CALC_DECO_CEILING */
+/*     Purpose: This subprogram calculates the deco ceiling (the safe ascent */
+/*     depth) in each compartment, based on the allowable gradients, and then */
+/*     finds the deepest deco ceiling across all compartments.  This deepest */
+/*     value (Deco Ceiling Depth) is then used by the Decompression Stop */
+/*     subroutine to determine the actual deco schedule. */
+/* =============================================================================== */
+
+int calc_deco_ceiling(float *deco_ceiling_depth,_Bool fallowable)
+{
+    /* System generated locals */
+    float r1, r2;
+    /* Local variables */
+    float weighted_allowable_gradient;
+    short i;
+    float compartment_deco_ceiling[16],
+    gas_loading,
+    tolerated_ambient_pressure;
+    float gradient_he, gradient_n2;
+
+    if(!vpm_b)
+        fallowable = true;
+    /* loop */
+    /* =============================================================================== */
+    /*     CALCULATIONS */
+    /*     Since there are two sets of allowable gradients being tracked, one for */
+    /*     helium and one for nitrogen, a "weighted allowable gradient" must be */
+    /*     computed each time based on the proportions of helium and nitrogen in */
+    /*     each compartment.  This proportioning follows the methodology of */
+    /*     Buhlmann/Keller.  If there is no helium and nitrogen in the compartment, */
+    /*     such as after extended periods of oxygen breathing, then the minimum value */
+    /*     across both gases will be used.  It is important to note that if a */
+    /*     compartment is empty of helium and nitrogen, then the weighted allowable */
+    /*     gradient formula cannot be used since it will result in division by zero. */
+    /* =============================================================================== */
+
+    for (i = 1; i <= 16; ++i)
+    {
+
+        // abfrage raus und pointer stattdessen
+        if(fallowable){
+            gradient_he = allowable_gradient_he[i-1];
+            gradient_n2 = allowable_gradient_n2[i-1];
+        }
+        else{
+            gradient_he = deco_gradient_he[i-1];
+            gradient_n2 = deco_gradient_n2[i-1];
+        }
+
+        gas_loading =	helium_pressure[i - 1] + nitrogen_pressure[i - 1];
+
+        if (gas_loading > 0)
+        {
+            weighted_allowable_gradient =
+            (gradient_he * helium_pressure[i - 1] +
+            gradient_n2 * nitrogen_pressure[i - 1]) /
+            (helium_pressure[i - 1] + nitrogen_pressure[i - 1]);
+
+            tolerated_ambient_pressure =
+            gas_loading +
+            CONSTANT_PRESSURE_OTHER_GASES -
+            weighted_allowable_gradient;
+        }
+        else
+        {
+            /* Computing MIN */
+            r1 = gradient_he;
+            r2 = gradient_n2;
+            weighted_allowable_gradient = fminf(r1,r2);
+
+            tolerated_ambient_pressure =
+            CONSTANT_PRESSURE_OTHER_GASES - weighted_allowable_gradient;
+        }
+
+        /* =============================================================================== */
+        /*     The tolerated ambient pressure cannot be less than zero absolute, i.e., */
+        /*     the vacuum of outer space! */
+        /* =============================================================================== */
+
+        if (tolerated_ambient_pressure < 0) {
+            tolerated_ambient_pressure = 0;
+        }
+        compartment_deco_ceiling[i - 1] =
+        tolerated_ambient_pressure - barometric_pressure;
+    }
+
+    /* =============================================================================== */
+    /*     The Deco Ceiling Depth is computed in a loop after all of the individual */
+    /*     compartment deco ceilings have been calculated.  It is important that the */
+    /*     Deco Ceiling Depth (max deco ceiling across all compartments) only be */
+    /*     extracted from the compartment values and not be compared against some */
+    /*     initialization value.  For example, if MAX(Deco_Ceiling_Depth . .) was */
+    /*     compared against zero, this could cause a program lockup because sometimes */
+    /*     the Deco Ceiling Depth needs to be negative (but not less than zero */
+    /*     absolute ambient pressure) in order to decompress to the last stop at zero */
+    /*     depth. */
+    /* =============================================================================== */
+
+    *deco_ceiling_depth = compartment_deco_ceiling[0];
+    for (i = 2; i <= 16; ++i)
+    {
+        /* Computing MAX */
+        r1 = *deco_ceiling_depth;
+        r2 = compartment_deco_ceiling[i - 1];
+        *deco_ceiling_depth = fmaxf(r1,r2);
+    }
+    return 0;
+} /* calc_deco_ceiling */
+
+
+
+/* =============================================================================== */
+/*     SUBROUTINE CALC_MAX_ACTUAL_GRADIENT */
+/*     Purpose: This subprogram calculates the actual supersaturation gradient */
+/*     obtained in each compartment as a result of the ascent profile during */
+/*     decompression.  Similar to the concept with crushing pressure, the */
+/*     supersaturation gradients are not cumulative over a multi-level, staged */
+/*     ascent.  Rather, it will be the maximum value obtained in any one discrete */
+/*     step of the overall ascent.  Thus, the program must compute and store the */
+/*     maximum actual gradient for each compartment that was obtained across all */
+/*     steps of the ascent profile.  This subroutine is invoked on the last pass */
+/*     through the deco stop loop block when the final deco schedule is being */
+/*     generated. */
+/*   */
+/*     The max actual gradients are later used by the VPM Repetitive Algorithm to */
+/*     determine if adjustments to the critical radii are required.  If the max */
+/*     actual gradient did not exceed the initial alllowable gradient, then no */
+/*     adjustment will be made.  However, if the max actual gradient did exceed */
+/*     the intitial allowable gradient, such as permitted by the Critical Volume */
+/*     Algorithm, then the critical radius will be adjusted (made larger) on the */
+/*     repetitive dive to compensate for the bubbling that was allowed on the */
+/*     previous dive.  The use of the max actual gradients is intended to prevent */
+/*     the repetitive algorithm from being overly conservative. */
+/* =============================================================================== */
+
+int calc_max_actual_gradient(float *deco_stop_depth)
+{
+    /* System generated locals */
+    float r1;
+
+    /* Local variables */
+    short i;
+    float compartment_gradient;
+
+    /* loop */
+    /* =============================================================================== */
+    /*     CALCULATIONS */
+    /*     Note: negative supersaturation gradients are meaningless for this */
+    /*     application, so the values must be equal to or greater than zero. */
+    /* =============================================================================== */
+
+    for (i = 1; i <= 16; ++i)
+    {
+        compartment_gradient =
+        helium_pressure[i - 1] +
+        nitrogen_pressure[i - 1] +
+        CONSTANT_PRESSURE_OTHER_GASES -
+        (*deco_stop_depth + barometric_pressure);
+        if (compartment_gradient <= 0.0f) {
+            compartment_gradient = 0.0f;
+        }
+        /* Computing MAX */
+        r1 = pVpm->max_actual_gradient[i - 1];
+        pVpm->max_actual_gradient[i - 1] = fmaxf(r1, compartment_gradient);
+    }
+    return 0;
+} /* calc_max_actual_gradient */
+
+/* =============================================================================== */
+/*     SUBROUTINE CALC_SURFACE_PHASE_VOLUME_TIME */
+/*     Purpose: This subprogram computes the surface portion of the total phase */
+/*     volume time.  This is the time factored out of the integration of */
+/*     supersaturation gradient x time over the surface interval.  The VPM */
+/*     considers the gradients that allow bubbles to form or to drive bubble */
+/*     growth both in the water and on the surface after the dive. */
+
+/*     This subroutine is a new development to the VPM algorithm in that it */
+/*     computes the time course of supersaturation gradients on the surface */
+/*     when both helium and nitrogen are present.  Refer to separate write-up */
+/*     for a more detailed explanation of this algorithm. */
+/* =============================================================================== */
+
+int calc_surface_phase_volume_time()
+{
+    /* Local variables */
+    float decay_time_to_zero_gradient;
+    short i;
+    float integral_gradient_x_time,
+    surface_inspired_n2_pressure;
+
+    /* loop */
+    /* =============================================================================== */
+    /*     CALCULATIONS */
+    /* =============================================================================== */
+
+    surface_inspired_n2_pressure =
+    (barometric_pressure - WATER_VAPOR_PRESSURE) * 0.79f;
+    for (i = 1; i <= 16; ++i)
+    {
+        if (nitrogen_pressure[i - 1] > surface_inspired_n2_pressure)
+        {
+            surface_phase_volume_time[i - 1] =
+            (helium_pressure[i - 1] / HELIUM_TIME_CONSTANT[i - 1] +
+            (nitrogen_pressure[i - 1] - surface_inspired_n2_pressure) /
+            NITROGEN_TIME_CONSTANT[i - 1]) /
+            (helium_pressure[i - 1] + nitrogen_pressure[i - 1] -
+            surface_inspired_n2_pressure);
+        } else if (nitrogen_pressure[i - 1] <= surface_inspired_n2_pressure &&
+        helium_pressure[i - 1] + nitrogen_pressure[i - 1] >= surface_inspired_n2_pressure)
+        {
+            decay_time_to_zero_gradient =
+            1.0f / (NITROGEN_TIME_CONSTANT[i - 1] - HELIUM_TIME_CONSTANT[i - 1]) *
+            log((surface_inspired_n2_pressure - nitrogen_pressure[i - 1]) /
+            helium_pressure[i - 1]);
+            integral_gradient_x_time =
+            helium_pressure[i - 1] /
+            HELIUM_TIME_CONSTANT[i - 1] *
+            (1.0f - expf(-HELIUM_TIME_CONSTANT[i - 1] *
+            decay_time_to_zero_gradient)) +
+            (nitrogen_pressure[i - 1] - surface_inspired_n2_pressure) /
+            NITROGEN_TIME_CONSTANT[i - 1] *
+            (1.0f - expf(-NITROGEN_TIME_CONSTANT[i - 1] *
+            decay_time_to_zero_gradient));
+            surface_phase_volume_time[i - 1] =
+            integral_gradient_x_time /
+            (helium_pressure[i - 1] +
+            nitrogen_pressure[i - 1] -
+            surface_inspired_n2_pressure);
+        } else {
+            surface_phase_volume_time[i - 1] = 0.0f;
+        }
+    }
+    return 0;
+} /* calc_surface_phase_volume_time */
+
+/* =============================================================================== */
+/*     SUBROUTINE CRITICAL_VOLUME */
+/*     Purpose: This subprogram applies the VPM Critical Volume Algorithm.  This */
+/*     algorithm will compute "relaxed" gradients for helium and nitrogen based */
+/*     on the setting of the Critical Volume Parameter Lambda. */
+/* =============================================================================== */
+
+int critical_volume(float *deco_phase_volume_time)
+{
+    /* System generated locals */
+    float r1;
+
+    /* Local variables */
+    float initial_allowable_grad_n2_pa,
+    initial_allowable_grad_he_pa,
+    parameter_lambda_pascals, b,
+    c;
+    short i;
+    float new_allowable_grad_n2_pascals,
+    phase_volume_time[16],
+    new_allowable_grad_he_pascals,
+    adj_crush_pressure_n2_pascals,
+    adj_crush_pressure_he_pascals;
+
+    /* loop */
+    /* =============================================================================== */
+    /*     CALCULATIONS */
+    /*     Note:  Since the Critical Volume Parameter Lambda was defined in units of */
+    /*     fsw-min in the original papers by Yount and colleauges, the same */
+    /*     convention is retained here.  Although Lambda is adjustable only in units */
+    /*     of fsw-min in the program settings (range from 6500 to 8300 with default */
+    /*     7500), it will convert to the proper value in Pascals-min in this */
+    /*     subroutine regardless of which diving pressure units are being used in */
+    /*     the main program - feet of seawater (fsw) or meters of seawater (msw). */
+    /*     The allowable gradient is computed using the quadratic formula (refer to */
+    /*     separate write-up posted on the Deco List web site). */
+    /* =============================================================================== */
+
+/**
+    ******************************************************************************
+    * @brief		critical_volume comment by hw
+    * @version 	V0.0.1
+    * @date   	19-April-2014
+    * @retval 	global: allowable_gradient_he[i], allowable_gradient_n2[i]
+    ******************************************************************************
+    */
+
+    parameter_lambda_pascals =
+    CRIT_VOLUME_PARAMETER_LAMBDA / 33.0f * 101325.0f;
+    for (i = 1; i <= 16; ++i)
+    {
+        phase_volume_time[i - 1] =
+        *deco_phase_volume_time + surface_phase_volume_time[i - 1];
+    }
+    for (i = 1; i <= 16; ++i)
+    {
+
+        adj_crush_pressure_he_pascals =
+        pVpm->adjusted_crushing_pressure_he[i - 1] / UNITS_FACTOR * 101325.0f;
+
+        initial_allowable_grad_he_pa =
+        pVpm->initial_allowable_gradient_he[i - 1] / UNITS_FACTOR * 101325.0f;
+
+        b = initial_allowable_grad_he_pa + parameter_lambda_pascals *
+        SURFACE_TENSION_GAMMA / (
+        SKIN_COMPRESSION_GAMMAC * phase_volume_time[i - 1]);
+
+        c = SURFACE_TENSION_GAMMA * (
+        SURFACE_TENSION_GAMMA * (
+        parameter_lambda_pascals * adj_crush_pressure_he_pascals)) /
+        (SKIN_COMPRESSION_GAMMAC *
+        (SKIN_COMPRESSION_GAMMAC * phase_volume_time[i - 1]));
+        /* Computing 2nd power */
+
+        r1 = b;
+
+        new_allowable_grad_he_pascals =
+        (b + sqrtf(r1 * r1 - c * 4.0f)) / 2.0f;
+
+        /* modify global variable */
+        allowable_gradient_he[i - 1] =
+        new_allowable_grad_he_pascals / 101325.0f * UNITS_FACTOR;
+    }
+
+    for (i = 1; i <= 16; ++i)
+    {
+        adj_crush_pressure_n2_pascals =
+        pVpm->adjusted_crushing_pressure_n2[i - 1] / UNITS_FACTOR * 101325.0f;
+
+        initial_allowable_grad_n2_pa =
+        pVpm->initial_allowable_gradient_n2[i - 1] / UNITS_FACTOR * 101325.0f;
+
+        b = initial_allowable_grad_n2_pa + parameter_lambda_pascals *
+        SURFACE_TENSION_GAMMA / (
+        SKIN_COMPRESSION_GAMMAC * phase_volume_time[i - 1]);
+
+        c = SURFACE_TENSION_GAMMA *
+        (SURFACE_TENSION_GAMMA *
+        (parameter_lambda_pascals * adj_crush_pressure_n2_pascals)) /
+        (SKIN_COMPRESSION_GAMMAC *
+        (SKIN_COMPRESSION_GAMMAC * phase_volume_time[i - 1]));
+        /* Computing 2nd power */
+
+        r1 = b;
+
+        new_allowable_grad_n2_pascals =
+        (b + sqrtf(r1 * r1 - c * 4.0f)) / 2.0f;
+
+        /* modify global variable */
+        allowable_gradient_n2[i - 1] =
+        new_allowable_grad_n2_pascals / 101325.0f * UNITS_FACTOR;
+    }
+    return 0;
+} /* critical_volume */
+
+/* =============================================================================== */
+/*     SUBROUTINE CALC_START_OF_DECO_ZONE */
+/*     Purpose: This subroutine uses the Bisection Method to find the depth at */
+/*     which the leading compartment just enters the decompression zone. */
+/*     Source:  "Numerical Recipes in Fortran 77", Cambridge University Press, */
+/*     1992. */
+/* =============================================================================== */
+
+int calc_start_of_deco_zone(float *starting_depth,
+                            float *rate,
+                            float *depth_start_of_deco_zone)
+{
+    /* Local variables */
+    float last_diff_change,
+    initial_helium_pressure,
+    mid_range_nitrogen_pressure;
+    short i, j;
+    float initial_inspired_n2_pressure,
+    cpt_depth_start_of_deco_zone,
+    low_bound,
+    initial_inspired_he_pressure,
+    high_bound_nitrogen_pressure,
+    nitrogen_rate,
+    function_at_mid_range,
+    function_at_low_bound,
+    high_bound,
+    mid_range_helium_pressure,
+    mid_range_time,
+    starting_ambient_pressure,
+    initial_nitrogen_pressure,
+    function_at_high_bound;
+
+    float time_to_start_of_deco_zone,
+    high_bound_helium_pressure,
+    helium_rate,
+    differential_change;
+    float fraction_helium_begin;
+    float fraction_helium_end;
+    float fraction_nitrogen_begin;
+    float fraction_nitrogen_end;
+    float ending_ambient_pressure;
+    float time_test;
+
+
+    /* loop */
+    /* =============================================================================== */
+    /*     CALCULATIONS */
+    /*     First initialize some variables */
+    /* =============================================================================== */
+
+    *depth_start_of_deco_zone = 0.0f;
+    starting_ambient_pressure =	*starting_depth + barometric_pressure;
+
+    //>>>>>>>>>>>>>>>>>>>>
+    //Test depth to calculate helium_rate and nitrogen_rate
+    ending_ambient_pressure = starting_ambient_pressure/2;
+
+    time_test = (ending_ambient_pressure - starting_ambient_pressure) / *rate;
+    decom_get_inert_gases(starting_ambient_pressure / 10, (&pDiveSettings->decogaslist[mix_number]), &fraction_nitrogen_begin, &fraction_helium_begin );
+    decom_get_inert_gases(ending_ambient_pressure   / 10, (&pDiveSettings->decogaslist[mix_number]), &fraction_nitrogen_end, &fraction_helium_end );
+    initial_inspired_he_pressure =	(starting_ambient_pressure - WATER_VAPOR_PRESSURE) * fraction_helium_begin;
+    initial_inspired_n2_pressure =	(starting_ambient_pressure - WATER_VAPOR_PRESSURE) * fraction_nitrogen_begin;
+    helium_rate = ((ending_ambient_pressure  - WATER_VAPOR_PRESSURE)* fraction_helium_end - initial_inspired_he_pressure)/time_test;
+    nitrogen_rate = ((ending_ambient_pressure  - WATER_VAPOR_PRESSURE)* fraction_nitrogen_end - initial_inspired_n2_pressure)/time_test;
+    //>>>>>>>>>>>>>>>>>>>>>
+    /*initial_inspired_he_pressure =
+    (starting_ambient_pressure - water_vapor_pressure) *
+    fraction_helium[mix_number - 1];
+    initial_inspired_n2_pressure =
+    (starting_ambient_pressure - water_vapor_pressure) *
+    fraction_nitrogen[mix_number - 1];
+    helium_rate = *rate * fraction_helium[mix_number - 1];
+    nitrogen_rate = *rate * fraction_nitrogen[mix_number - 1];*/
+
+    /* =============================================================================== */
+    /*     ESTABLISH THE BOUNDS FOR THE ROOT SEARCH USING THE BISECTION METHOD */
+    /*     AND CHECK TO MAKE SURE THAT THE ROOT WILL BE WITHIN BOUNDS.  PROCESS */
+    /*     EACH COMPARTMENT INDIVIDUALLY AND FIND THE MAXIMUM DEPTH ACROSS ALL */
+    /*     COMPARTMENTS (LEADING COMPARTMENT) */
+    /*     In this case, we are solving for time - the time when the gas tension in */
+    /*     the compartment will be equal to ambient pressure.  The low bound for time */
+    /*     is set at zero and the high bound is set at the time it would take to */
+    /*     ascend to zero ambient pressure (absolute).  Since the ascent rate is */
+    /*     negative, a multiplier of -1.0 is used to make the time positive.  The */
+    /*     desired point when gas tension equals ambient pressure is found at a time */
+    /*     somewhere between these endpoints.  The algorithm checks to make sure that */
+    /*     the solution lies in between these bounds by first computing the low bound */
+    /*     and high bound function values. */
+    /* =============================================================================== */
+
+    low_bound = 0.;
+    high_bound = starting_ambient_pressure / *rate * -1.0f;
+    for (i = 1; i <= 16; ++i)
+    {
+        initial_helium_pressure = helium_pressure[i - 1];
+        initial_nitrogen_pressure = nitrogen_pressure[i - 1];
+        function_at_low_bound =
+        initial_helium_pressure +
+        initial_nitrogen_pressure +
+        CONSTANT_PRESSURE_OTHER_GASES -
+        starting_ambient_pressure;
+        high_bound_helium_pressure =
+        schreiner_equation__2(&initial_inspired_he_pressure,
+        &helium_rate,
+        &high_bound,
+        &HELIUM_TIME_CONSTANT[i - 1],
+        &initial_helium_pressure);
+        high_bound_nitrogen_pressure =
+        schreiner_equation__2(&initial_inspired_n2_pressure,
+        &nitrogen_rate,
+        &high_bound,
+        &NITROGEN_TIME_CONSTANT[i - 1],
+        &initial_nitrogen_pressure);
+        function_at_high_bound = high_bound_helium_pressure +
+        high_bound_nitrogen_pressure +
+        CONSTANT_PRESSURE_OTHER_GASES;
+        if (function_at_high_bound * function_at_low_bound >= 0.0f)
+        {
+            printf("\nERROR! ROOT IS NOT WITHIN BRACKETS");
+        }
+
+        /* =============================================================================== */
+        /*     APPLY THE BISECTION METHOD IN SEVERAL ITERATIONS UNTIL A SOLUTION WITH */
+        /*     THE DESIRED ACCURACY IS FOUND */
+        /*     Note: the program allows for up to 100 iterations.  Normally an exit will */
+        /*     be made from the loop well before that number.  If, for some reason, the */
+        /*     program exceeds 100 iterations, there will be a pause to alert the user. */
+        /* =============================================================================== */
+
+        if (function_at_low_bound < 0.0f)
+        {
+            time_to_start_of_deco_zone = low_bound;
+            differential_change = high_bound - low_bound;
+        } else {
+            time_to_start_of_deco_zone = high_bound;
+            differential_change = low_bound - high_bound;
+        }
+        for (j = 1; j <= 100; ++j)
+        {
+            last_diff_change = differential_change;
+            differential_change = last_diff_change * 0.5f;
+            mid_range_time =
+            time_to_start_of_deco_zone +
+            differential_change;
+            mid_range_helium_pressure =
+            schreiner_equation__2(&initial_inspired_he_pressure,
+            &helium_rate,
+            &mid_range_time,
+            &HELIUM_TIME_CONSTANT[i - 1],
+            &initial_helium_pressure);
+            mid_range_nitrogen_pressure =
+            schreiner_equation__2(&initial_inspired_n2_pressure,
+            &nitrogen_rate,
+            &mid_range_time,
+            &NITROGEN_TIME_CONSTANT[i - 1],
+            &initial_nitrogen_pressure);
+            function_at_mid_range =
+            mid_range_helium_pressure +
+            mid_range_nitrogen_pressure +
+            CONSTANT_PRESSURE_OTHER_GASES -
+            (starting_ambient_pressure + *rate * mid_range_time);
+            if (function_at_mid_range <= 0.0f) {
+                time_to_start_of_deco_zone = mid_range_time;
+            }
+            if( fabs(differential_change) < 0.001f
+             || function_at_mid_range == 0.0f)
+            {
+                goto L170;
+            }
+            /* L150: */
+        }
+        printf("\nERROR! ROOT SEARCH EXCEEDED MAXIMUM ITERATIONS");
+        //pause();
+
+        /* =============================================================================== */
+        /*     When a solution with the desired accuracy is found, the program jumps out */
+        /*     of the loop to Line 170 and assigns the solution value for the individual */
+        /*     compartment. */
+        /* =============================================================================== */
+
+        L170:
+        cpt_depth_start_of_deco_zone =
+        starting_ambient_pressure +
+        *rate * time_to_start_of_deco_zone -
+        barometric_pressure;
+
+        /* =============================================================================== */
+        /*     The overall solution will be the compartment with the maximum depth where */
+        /*     gas tension equals ambient pressure (leading compartment). */
+        /* =============================================================================== */
+
+        *depth_start_of_deco_zone =
+        fmaxf(*depth_start_of_deco_zone, cpt_depth_start_of_deco_zone);
+        /* L200: */
+    }
+    return 0;
+} /* calc_start_of_deco_zone */
+
+/* =============================================================================== */
+/*     SUBROUTINE PROJECTED_ASCENT */
+/*     Purpose: This subprogram performs a simulated ascent outside of the main */
+/*     program to ensure that a deco ceiling will not be violated due to unusual */
+/*     gas loading during ascent (on-gassing).  If the deco ceiling is violated, */
+/*     the stop depth will be adjusted deeper by the step size until a safe */
+/*     ascent can be made. */
+/* =============================================================================== */
+
+int projected_ascent(float *starting_depth,
+                     float *rate,
+                     float *deco_stop_depth,
+                     float *step_size)
+{
+    /* Local variables */
+    float weighted_allowable_gradient,
+    ending_ambient_pressure,
+    temp_gas_loading[16];
+    int i;
+    float allowable_gas_loading[16];
+    float temp_nitrogen_pressure[16];
+    float temp_helium_pressure[16];
+    float run_time_save = 0;
+
+    /* loop */
+    /* =============================================================================== */
+    /*     CALCULATIONS */
+    /* =============================================================================== */
+
+
+    L665:
+    ending_ambient_pressure = *deco_stop_depth + barometric_pressure;
+    for (i = 1; i <= 16; ++i) {
+        temp_helium_pressure[i - 1] = helium_pressure[i - 1];
+        temp_nitrogen_pressure[i - 1] = nitrogen_pressure[i - 1];
+    }
+    run_time_save = run_time;
+    gas_loadings_ascent_descen(temp_helium_pressure, temp_nitrogen_pressure, *starting_depth,*deco_stop_depth,*rate,true);
+    run_time = run_time_save;
+
+    for (i = 1; i <= 16; ++i)
+    {
+        temp_gas_loading[i - 1] =
+        temp_helium_pressure[i - 1] +
+        temp_nitrogen_pressure[i - 1];
+        if (temp_gas_loading[i - 1] > 0.0f)
+        {
+            weighted_allowable_gradient =
+            (allowable_gradient_he[i - 1] *
+            temp_helium_pressure[i - 1] +
+            allowable_gradient_n2[i - 1] *
+            temp_nitrogen_pressure[i - 1]) / temp_gas_loading[i - 1];
+        } else {
+            /* Computing MIN */
+            weighted_allowable_gradient = fminf(allowable_gradient_he[i - 1],allowable_gradient_n2[i - 1]);
+        }
+        allowable_gas_loading[i - 1] =
+        ending_ambient_pressure +
+        weighted_allowable_gradient -
+        CONSTANT_PRESSURE_OTHER_GASES;
+        /* L670: */
+    }
+    for (i = 1; i <= 16; ++i) {
+        if (temp_gas_loading[i - 1] > allowable_gas_loading[i - 1]) {
+            *deco_stop_depth += *step_size;
+            goto L665;
+        }
+        /* L671: */
+    }
+    return 0;
+} /* projected_ascent */
+
+/* =============================================================================== */
+/*     SUBROUTINE DECOMPRESSION_STOP */
+/*     Purpose: This subprogram calculates the required time at each */
+/*     decompression stop. */
+/* =============================================================================== */
+
+void decompression_stop(float *deco_stop_depth,
+                        float *step_size,
+                        _Bool final_deco_calculation)
+{
+    /* Local variables */
+    float inspired_nitrogen_pressure;
+//	short last_segment_number;
+//	float weighted_allowable_gradient;
+    float initial_helium_pressure[16];
+ /* by hw */
+    float initial_CNS;
+
+    //static float time_counter;
+    short i;
+    float ambient_pressure;
+    float inspired_helium_pressure,
+    next_stop;
+    //last_run_time,
+    //temp_segment_time;
+
+    float deco_ceiling_depth,
+    initial_nitrogen_pressure[16];
+    //round_up_operation;
+    float fraction_helium_begin;
+    float fraction_nitrogen_begin;
+    int count = 0;
+    _Bool buehlmann_wait = false;
+    float tissue_He_saturation[16];
+    float tissue_N2_saturation[16];
+    float vpm_buehlmann_safety_gradient = 1.0f - (((float)pDiveSettings->vpm_conservatism) / 40);
+    /* loop */
+    /* =============================================================================== */
+    /*     CALCULATIONS */
+    /* =============================================================================== */
+
+    segment_time = 0;
+//	temp_segment_time = segment_time;
+    ambient_pressure = *deco_stop_depth + barometric_pressure;
+    //ending_ambient_pressure = ambient_pressure;
+    decom_get_inert_gases(ambient_pressure / 10, (&pDiveSettings->decogaslist[mix_number]), &fraction_nitrogen_begin, &fraction_helium_begin );
+
+    if(*deco_stop_depth == (float)(pDiveSettings->last_stop_depth_bar * 10))
+        next_stop = 0;
+    else
+    {
+        next_stop = *deco_stop_depth - *step_size;
+        next_stop = fmaxf(next_stop,(float)pDiveSettings->last_stop_depth_bar * 10);
+    }
+
+    inspired_helium_pressure =
+    (ambient_pressure - WATER_VAPOR_PRESSURE) * fraction_helium_begin;
+    inspired_nitrogen_pressure =
+    (ambient_pressure - WATER_VAPOR_PRESSURE) * fraction_nitrogen_begin;
+
+    /* =============================================================================== */
+    /*     Check to make sure that program won't lock up if unable to decompress */
+    /*     to the next stop.  If so, write error message and terminate program. */
+    /* =============================================================================== */
+
+    //deco_ceiling_depth = next_stop +1;	//deco_ceiling_depth = next_stop + 1;
+    if(!vpm_violates_buehlmann)
+        calc_deco_ceiling(&deco_ceiling_depth, false); //weg, weil auf jeden Fall schleife für safety und so konservativer
+    else
+        deco_ceiling_depth = next_stop + 1;
+
+        if(deco_ceiling_depth > next_stop)
+        {
+            while (deco_ceiling_depth > next_stop)
+            {
+
+                segment_time  +=  60;
+                if(segment_time >= 999 )
+                {
+                    segment_time = 999 ;
+                    run_time += segment_time;
+                    return;
+                }
+                //goto L700;
+                initial_CNS = gCNS_VPM;
+                decom_oxygen_calculate_cns_exposure(60*60,&pDiveSettings->decogaslist[mix_number],ambient_pressure/10,&gCNS_VPM);
+                for (i = 0; i < 16; i++)
+                {
+                    initial_helium_pressure[i] = helium_pressure[i];
+                    initial_nitrogen_pressure[i] = nitrogen_pressure[i];
+                    helium_pressure[i] 		+=	(inspired_helium_pressure - helium_pressure[i]) 		*	float_buehlmann_He_factor_expositon_one_hour[i];
+                    nitrogen_pressure[i]	+= 	(inspired_nitrogen_pressure - nitrogen_pressure[i]) *	float_buehlmann_N2_factor_expositon_one_hour[i];
+                }
+                calc_deco_ceiling(&deco_ceiling_depth, false);
+            }
+            if(deco_ceiling_depth < next_stop)
+            {
+                segment_time  -=  60;
+                gCNS_VPM = initial_CNS;
+                for (i = 0; i < 16; i++)
+                {
+                    helium_pressure[i] = initial_helium_pressure[i];
+                    nitrogen_pressure[i] = initial_nitrogen_pressure[i];
+                }
+                deco_ceiling_depth = next_stop +1;
+            }
+            count = 0;
+            while (deco_ceiling_depth > next_stop && count < 13)
+            {
+                count++;
+                segment_time  +=  5;
+                //goto L700;
+                initial_CNS = gCNS_VPM;
+                decom_oxygen_calculate_cns_exposure(60*5,&pDiveSettings->decogaslist[mix_number],ambient_pressure/10,&gCNS_VPM);
+                for (i = 0; i < 16; i++)
+                {
+                    initial_helium_pressure[i] = helium_pressure[i];
+                    initial_nitrogen_pressure[i] = nitrogen_pressure[i];
+                    helium_pressure[i] 		+=	(inspired_helium_pressure - helium_pressure[i]) 		*	float_buehlmann_He_factor_expositon_five_minutes[i];
+                    nitrogen_pressure[i]	+= 	(inspired_nitrogen_pressure - nitrogen_pressure[i]) *	float_buehlmann_N2_factor_expositon_five_minutes[i];
+                }
+                calc_deco_ceiling(&deco_ceiling_depth, false);
+            }
+            if(deco_ceiling_depth < next_stop)
+            {
+                segment_time  -=  5;
+                gCNS_VPM = initial_CNS;
+                for (i = 0; i < 16; i++) {
+                    helium_pressure[i] = initial_helium_pressure[i];
+                    nitrogen_pressure[i] = initial_nitrogen_pressure[i];
+                }
+                deco_ceiling_depth = next_stop +1;
+            }
+            buehlmann_wait = false;
+            while (buehlmann_wait || (deco_ceiling_depth > next_stop))
+            {
+                //time_counter = temp_segment_time;
+                segment_time  +=  1;
+
+                if(segment_time >= 999 )
+                {
+                    segment_time = 999 ;
+                    run_time += segment_time;
+                    return;
+                }
+                //goto L700;
+                initial_CNS = gCNS_VPM;
+                decom_oxygen_calculate_cns_exposure(60*1,&pDiveSettings->decogaslist[mix_number],ambient_pressure/10,&gCNS_VPM);
+                for (i = 0; i < 16; i++)
+                {
+                    initial_helium_pressure[i] = helium_pressure[i];
+                    initial_nitrogen_pressure[i] = nitrogen_pressure[i];
+                    helium_pressure[i] 		+=	(inspired_helium_pressure - helium_pressure[i]) 		*	float_buehlmann_He_factor_expositon_one_minute[i];
+                    nitrogen_pressure[i]	+= 	(inspired_nitrogen_pressure - nitrogen_pressure[i]) *	float_buehlmann_N2_factor_expositon_one_minute[i];
+                }
+                if(!buehlmann_wait)
+                    calc_deco_ceiling(&deco_ceiling_depth, false);
+
+                if(buehlmannSafety && final_deco_calculation && !(deco_ceiling_depth > next_stop))
+                {
+                    for (i = 0; i < 16; i++)
+                    {
+                        tissue_He_saturation[i] = helium_pressure[i] / 10;
+                        tissue_N2_saturation[i] = nitrogen_pressure[i] / 10;
+                    }
+                    if( (fabsf(nitrogen_pressure[15] - inspired_nitrogen_pressure) < 0.00001f) && (fabsf(helium_pressure[15] - inspired_helium_pressure) < 0.00001f)
+                     && (fabsf(nitrogen_pressure[0] - inspired_nitrogen_pressure) < 0.00001f) && (fabsf(helium_pressure[0] - inspired_helium_pressure) < 0.00001f))
+                    {
+                         buehlmann_wait_exceeded = true;
+                        break;
+                    }
+
+                    if(decom_tissue_test_tolerance(tissue_N2_saturation, tissue_He_saturation, vpm_buehlmann_safety_gradient, (next_stop / 10.0f) + pInput->pressure_surface_bar))
+                        break;
+
+                    buehlmann_wait = true;
+                }
+            }
+            if(buehlmann_wait)
+                vpm_violates_buehlmann = true;
+
+                if(!buehlmann_wait)
+            {
+                if(deco_ceiling_depth < next_stop)
+                {
+                    segment_time  -=  1;
+                    gCNS_VPM = initial_CNS;
+                    for (i = 0; i < 16; i++) {
+                        helium_pressure[i] = initial_helium_pressure[i];
+                        nitrogen_pressure[i] = initial_nitrogen_pressure[i];
+                    }
+                    deco_ceiling_depth = next_stop +1;
+                }
+                while (deco_ceiling_depth > next_stop)
+                {
+                    //time_counter = temp_segment_time;
+                    segment_time  +=  (float) 1.0f / 3.0f;
+                    //goto L700;
+                    initial_CNS = gCNS_VPM;
+                    decom_oxygen_calculate_cns_exposure(20,&pDiveSettings->decogaslist[mix_number],ambient_pressure/10,&gCNS_VPM);
+                    for (i = 0; i < 16; i++)
+                    {
+                        helium_pressure[i] 		+=	(inspired_helium_pressure - helium_pressure[i]) 		*	float_buehlmann_He_factor_expositon_20_seconds[i];
+                        nitrogen_pressure[i]	+= 	(inspired_nitrogen_pressure - nitrogen_pressure[i]) *	float_buehlmann_N2_factor_expositon_20_seconds[i];
+                    }
+                    calc_deco_ceiling(&deco_ceiling_depth, false);
+                }
+            }
+    }
+
+/*float pressure_save =dive_data.pressure;
+dive_data.pressure = ambient_pressure/10;
+tissues_exposure_stage(st_deco_test,(int)(segment_time * 60), &dive_data, &gaslist);
+dive_data.pressure = pressure_save;*/
+    run_time += segment_time;
+    return;
+} /* decompression_stop */
+
+/* =============================================================================== */
+// SUROUTINE BOYLES_LAW_COMPENSATION
+//     Purpose: This subprogram calculates the reduction in allowable gradients
+//    with decreasing ambient pressure during the decompression profile based
+//    on Boyle's Law considerations.
+//===============================================================================
+void BOYLES_LAW_COMPENSATION (float* First_Stop_Depth,
+                              float*  Deco_Stop_Depth,
+                              float* Step_Size)
+{
+    short i;
+
+    float Next_Stop;
+    float Ambient_Pressure_First_Stop, Ambient_Pressure_Next_Stop;
+    float Amb_Press_First_Stop_Pascals, Amb_Press_Next_Stop_Pascals;
+    float A, B, C, Low_Bound, High_Bound, Ending_Radius;
+    float Deco_Gradient_Pascals;
+    float Allow_Grad_First_Stop_He_Pa, Radius_First_Stop_He;
+    float Allow_Grad_First_Stop_N2_Pa, Radius_First_Stop_N2;
+
+    //===============================================================================
+    //     LO//AL ARRAYS
+    //===============================================================================
+//	float Radius1_He[16], Radius2_He[16];
+//	float Radius1_N2[16], Radius2_N2[16];
+    float root_factor;
+
+    //===============================================================================
+    //     CALCULATIONS
+    //===============================================================================
+    Next_Stop = *Deco_Stop_Depth - *Step_Size;
+
+    Ambient_Pressure_First_Stop = *First_Stop_Depth +
+                                 barometric_pressure;
+
+    Ambient_Pressure_Next_Stop = Next_Stop + barometric_pressure;
+
+    Amb_Press_First_Stop_Pascals =  (Ambient_Pressure_First_Stop/UNITS_FACTOR) * 101325.0f;
+
+    Amb_Press_Next_Stop_Pascals =
+            (Ambient_Pressure_Next_Stop/UNITS_FACTOR) * 101325.0f;
+    root_factor = powf(Amb_Press_First_Stop_Pascals/Amb_Press_Next_Stop_Pascals,1.0f / 3.0f);
+
+    for( i = 0; i < 16;i++)
+    {
+            Allow_Grad_First_Stop_He_Pa =
+                        (allowable_gradient_he[i]/UNITS_FACTOR) * 101325.0f;
+
+            Radius_First_Stop_He = (2.0f * SURFACE_TENSION_GAMMA) /
+                                    Allow_Grad_First_Stop_He_Pa;
+
+//			Radius1_He[i] = Radius_First_Stop_He;
+            A = Amb_Press_Next_Stop_Pascals;
+            B = -2.0f * SURFACE_TENSION_GAMMA;
+            C = (Amb_Press_First_Stop_Pascals + (2.0f * SURFACE_TENSION_GAMMA)/
+                     Radius_First_Stop_He)* Radius_First_Stop_He*
+                     (Radius_First_Stop_He*(Radius_First_Stop_He));
+            Low_Bound = Radius_First_Stop_He;
+            High_Bound = Radius_First_Stop_He * root_factor;
+                        //*pow(Amb_Press_First_Stop_Pascals/Amb_Press_Next_Stop_Pascals,1.0/3.0);
+                    //*(Amb_Press_First_Stop_Pascals/Amb_Press_Next_Stop_Pascals)**(1.0/3.0);
+
+            radius_root_finder(&A,&B,&C, &Low_Bound, &High_Bound,
+                                             &Ending_Radius);
+
+//			Radius2_He[i] = Ending_Radius;
+            Deco_Gradient_Pascals = (2.0f * SURFACE_TENSION_GAMMA) /
+                                    Ending_Radius;
+
+            deco_gradient_he[i] = (Deco_Gradient_Pascals / 101325.0f)*
+                                     UNITS_FACTOR;
+
+    }
+
+    for( i = 0; i < 16;i++)
+    {
+        Allow_Grad_First_Stop_N2_Pa =
+                             (allowable_gradient_n2[i]/UNITS_FACTOR) * 101325.0f;
+
+        Radius_First_Stop_N2 = (2.0f * SURFACE_TENSION_GAMMA) /
+                                                         Allow_Grad_First_Stop_N2_Pa;
+
+//		Radius1_N2[i] = Radius_First_Stop_N2;
+        A = Amb_Press_Next_Stop_Pascals;
+        B = -2.0f * SURFACE_TENSION_GAMMA;
+        C = (Amb_Press_First_Stop_Pascals + (2.0f * SURFACE_TENSION_GAMMA)/
+                        Radius_First_Stop_N2)* Radius_First_Stop_N2*
+                        (Radius_First_Stop_N2*(Radius_First_Stop_N2));
+        Low_Bound = Radius_First_Stop_N2;
+        High_Bound = Radius_First_Stop_N2* root_factor;//pow(Amb_Press_First_Stop_Pascals/Amb_Press_Next_Stop_Pascals,1.0/3.0);
+
+        //High_Bound = Radius_First_Stop_N2*exp(log(Amb_Press_First_Stop_Pascals/Amb_Press_Next_Stop_Pascals)/3);
+        radius_root_finder(&A,&B,&C, &Low_Bound, &High_Bound,
+                                                                        &Ending_Radius);
+
+//		Radius2_N2[i] = Ending_Radius;
+        Deco_Gradient_Pascals = (2.0f * SURFACE_TENSION_GAMMA) /
+                                                         Ending_Radius;
+
+        deco_gradient_n2[i] = (Deco_Gradient_Pascals / 101325.0f)*
+                                                        UNITS_FACTOR;
+    }
+}
+
+/* =============================================================================== */
+//	vpm_calc_nullzeit
+//     Purpose: This function calcs zero time (time where no decostops are needed)
+//===============================================================================
+int vpm_calc_nullzeit(void)
+{
+    static float future_helium_pressure[16];
+    static float future_nitrogen_pressure[16];
+    static int temp_segment_time;
+    static int mix_number;
+    static float inspired_helium_pressure;
+    static float inspired_nitrogen_pressure;
+
+    float previous_helium_pressure[16];
+    float previous_nitrogen_pressure[16];
+    float ambient_pressure;
+    float fraction_helium_begin;
+    float fraction_nitrogen_begin;
+    int i = 0;
+    int count = 0;
+    int status = CALC_END;
+    //if(begin)
+    //{
+        for(i = 0; i < 16;i++)
+        {
+            future_helium_pressure[i] =  pInput->tissue_helium_bar[i] * 10;//tissue_He_saturation[st_dive][i] * 10;
+            future_nitrogen_pressure[i] = pInput->tissue_nitrogen_bar[i] * 10;
+        }
+        temp_segment_time = 0;
+
+        mix_number = 0;
+        ambient_pressure = pInput->pressure_ambient_bar  * 10;
+//		 fraction_helium_begin;
+//		 fraction_nitrogen_begin;
+        decom_get_inert_gases( ambient_pressure / 10,  (&pDiveSettings->decogaslist[mix_number]) , &fraction_nitrogen_begin, &fraction_helium_begin );
+        inspired_helium_pressure =(ambient_pressure - WATER_VAPOR_PRESSURE) * fraction_helium_begin;
+        inspired_nitrogen_pressure =(ambient_pressure - WATER_VAPOR_PRESSURE) *fraction_nitrogen_begin;
+
+        //if(!nullzeit_unter60)
+        //{
+            status = CALC_END;
+            while (status == CALC_END)
+            {
+                count++;
+                //if(count == 7)
+                    //return CALC_NULLZEIT2;
+                temp_segment_time  +=  60;
+                if(temp_segment_time >= 300)
+                {
+                    pDecoInfo->output_ndl_seconds = temp_segment_time * 60;
+                    return CALC_NULLZEIT;
+                }
+                run_time += 60;
+                //goto L700;
+                for (i = 1; i <= 16; ++i) {
+                    previous_helium_pressure[i-1] = future_helium_pressure[i - 1];
+                    previous_nitrogen_pressure[i - 1] = future_nitrogen_pressure[i - 1];
+                    future_helium_pressure[i - 1] = future_helium_pressure[i - 1] + (inspired_helium_pressure - future_helium_pressure[i - 1]) * float_buehlmann_He_factor_expositon_one_hour[i-1];
+                    future_nitrogen_pressure[i - 1] = future_nitrogen_pressure[i - 1] + (inspired_nitrogen_pressure - future_nitrogen_pressure[i - 1]) * float_buehlmann_N2_factor_expositon_one_hour[i-1];
+                    helium_pressure[i - 1] = future_helium_pressure[i - 1];
+                    nitrogen_pressure[i - 1] = future_nitrogen_pressure[i - 1];
+                }
+                vpm_calc_deco();
+                while((status = vpm_calc_critcal_volume(true,true)) == CALC_CRITICAL);
+
+            }
+
+            temp_segment_time  -=  60;
+            run_time -= 60;
+            for (i = 1; i <= 16; ++i)
+            {
+                future_helium_pressure[i - 1] = previous_helium_pressure[i-1];
+                future_nitrogen_pressure[i - 1] = previous_nitrogen_pressure[i - 1];
+            }
+        //}
+    //}
+    //if(!nullzeit_unter60 || begin ||  temp_segment_time > 10)
+    //{
+
+        status = CALC_END;
+        if(temp_segment_time < 60)
+            nullzeit_unter60 = true;
+
+        while (status == CALC_END)
+        {
+        //	count++;
+            //if(count >= 5)
+                //return CALC_NULLZEIT2;
+            temp_segment_time  +=  5;
+            if(temp_segment_time >= 300)
+            {
+                pDecoInfo->output_ndl_seconds = temp_segment_time * 60;
+                return CALC_NULLZEIT;
+            }
+            if(nullzeit_unter60 && temp_segment_time > 60)
+            {
+                nullzeit_unter60 = false;
+                //tts[NULLZEIT] = temp_segment_time * 60;
+                return CALC_NULLZEIT;
+            }
+            run_time += 5;
+            //goto L700;
+            for (i = 1; i <= 16; ++i) {
+                previous_helium_pressure[i-1] = future_helium_pressure[i - 1];
+                previous_nitrogen_pressure[i - 1] = future_nitrogen_pressure[i - 1];
+                future_helium_pressure[i - 1] = future_helium_pressure[i - 1] + (inspired_helium_pressure - future_helium_pressure[i - 1]) * float_buehlmann_He_factor_expositon_five_minutes[i-1];
+                future_nitrogen_pressure[i - 1] = future_nitrogen_pressure[i - 1] + (inspired_nitrogen_pressure - future_nitrogen_pressure[i - 1]) * float_buehlmann_N2_factor_expositon_five_minutes[i-1];
+                helium_pressure[i - 1] = future_helium_pressure[i - 1];
+                nitrogen_pressure[i - 1] = future_nitrogen_pressure[i - 1];
+            }
+            vpm_calc_deco();
+            while((status =vpm_calc_critcal_volume(true,true)) == CALC_CRITICAL);
+        }
+        temp_segment_time  -=  5;
+        run_time -= 5;
+        for (i = 1; i <= 16; ++i) {
+            future_helium_pressure[i - 1] = previous_helium_pressure[i-1];
+            future_nitrogen_pressure[i - 1] = previous_nitrogen_pressure[i - 1];
+        }
+        status = CALC_END;
+        //if(temp_segment_time < 5)
+            //count = 2;
+    //}
+    //else
+        //count = 1;
+    if(temp_segment_time <= 20)
+    {
+        while (status == CALC_END)
+        {
+            //time_counter = temp_segment_time;
+            //count++;
+            //if(count > 2)
+                //return CALC_NULLZEIT2;
+            temp_segment_time  +=  minimum_deco_stop_time;
+            run_time += minimum_deco_stop_time;
+            //goto L700;
+            for (i = 1; i <= 16; ++i) {
+                future_helium_pressure[i - 1] = future_helium_pressure[i - 1] + (inspired_helium_pressure - future_helium_pressure[i - 1]) *	float_buehlmann_He_factor_expositon_one_minute[i-1];
+                future_nitrogen_pressure[i - 1] = future_nitrogen_pressure[i - 1] + (inspired_nitrogen_pressure - future_nitrogen_pressure[i - 1]) * float_buehlmann_N2_factor_expositon_one_minute[i-1];
+                helium_pressure[i - 1] = future_helium_pressure[i - 1];
+                nitrogen_pressure[i - 1] =future_nitrogen_pressure[i - 1];
+
+            }
+            vpm_calc_deco();
+            while((status =vpm_calc_critcal_volume(true,true)) == CALC_CRITICAL);
+
+        }
+    }
+    else
+        temp_segment_time += 5;
+    pDecoInfo->output_ndl_seconds = temp_segment_time * 60;
+    if(temp_segment_time > 1)
+        return CALC_NULLZEIT;
+    else
+        return CALC_BEGIN;
+}