firelib.c

这是一个GPS相关的程序

/*
 *******************************************************************************
 *
 *  fireLib.c
 *
 *  Description
 *      Library of BEHAVE (Andrews 1986) fire behavior algorithms
 *      encapsulated and optimized for fire behavior simulation.
 *
 *  Legalities
 *      Copyright (c) 1996 Collin D. Bevins.
 *      See the file "license.txt" for information on usage and
 *      redistribution of this file, and for a DISCLAIMER OF ALL WARRANTIES.
 *
 *  Naming Conventions
 *      All function names begin with "Fire_".
 *      All fuel model and behavior parameter access macros begin with "Fuel_".
 *      All fuel catalog parameter access macros begin with "FuelCat_".
 *
 *  Functions
 *      There are 8 functions to create and destroy fuel models and catalogs:
 *
 *          Fire_FuelCatalogCreate(name, maxModels)
 *              Creates a new fuel catalog capable of holding maxModels.
 *
 *          Fire_FuelCatalogCreateStandard(name, maxModels)
 *              Creates a new fuel catalog capable of holding maxModels,
 *              and fills models 0-13 with standard fire behavior models.
 *
 *          Fire_FuelModelCreate(catalog, model, name, desc, depth, mext,
 *                  adjust, maxParticles)
 *              Adds or replaces a fuel model in the catalog.  The model will
 *              accept up to maxParticles particles.
 *
 *          Fire_FuelModelExists(catalog, model)
 *              Returns 1 if model exists within the catalog.
 *
 *          Fire_FuelParticleAdd(catalog, model, live, load, savr, dens, heat,
 *                  stot, seff)
 *              Adds a fuel particle to a fuel model.
 *
 *          Fire_FlameLengthTable ( catalog, flameClasses, flameStep )
 *              Creates a flame length look-up table containing flameClasses
 *              number of classes, with each class spanning "flameStep"
 *              feet of flame length.  Creating a flame length table can
 *              significantly improve performance at the expense of user
 *              specified precision.
 *
 *          Fire_FuelModelDestroy(catalog, model)
 *              Destroys the model within the catalog.
 *
 *          Fire_FuelCatalogDestroy(catalog)
 *              Destroys the catalog and all models within it.
 *
 *      There are 5 functions to process data within fuel models:
 *
 *          Fire_FuelCombustion(catalog, model)
 *              Computes all the fuel-dependent model variables.
 *              Called only once for each fuel model.
 *              Called automatically by Fire_SpreadNoWindNoSlope().
 *
 *          Fire_SpreadNoWindNoSlope(catalog, model, moisture[6])
 *              Determines reaction intensity, heat per unit area, and the
 *              no-wind no-slope spread rate.
 *
 *          Fire_SpreadWindSlopeMax(catalog, model, windFpm, windDeg, slope,
 *                  aspectDeg)
 *              Determines maximum spread rate and azimuth of maximum spread
 *              based upon input parameters and results of the most recent
 *              call to Fire_SpreadNoWindNoSlope() for this model.
 *
 *          Fire_SpreadAtAzimuth(catalog, model, azimuth, doWhich)
 *              Determines the spread rate in the specified azimuth based
 *              upon the results of the most recent call to
 *              Fire_SpreadWindSlopeMax() for this model.  The "doWhich"
 *              parameter is the result of ORing the constants FIRE_BYRAMS,
 *              FIRE_FLAME, and FIRE_SCORCH to request computation of the
 *              associated fire variables.
 *
 *          Fire_FlameScorch(catalog, model, doWhich)
 *              Determines the flame length and/or scorch height based upon
 *              the most recent call to Fire_SpreadAtAzimuth().
 *
 *  History
 *      1996/09/04  Version 1.0.0 release.
 *      1999/03/05  Fixed NNFL07 live SAVR from 1500 to 1550.
 *
 *******************************************************************************
 */

#include "fireLib.h"

#ifndef M_PI
#define M_PI 3.14159
#endif

/*
 *******************************************************************************
 *
 *  Fire_FuelCombustion()
 *
 *  Description
 *      Calculates and stores all the fuel-dependent combustion variables.
 *
 *  Side Effects
 *      All combustion varaiables are reclaculated for the model.
 *      All behavior and environment variables are reset to zero.
 *
 *  Function Returns
 *      FIRE_STATUS_OK or FIRE_STATUS_ERROR.
 *      Return status and error text are stored in the Fire Catalog's buffers.
 *
 *******************************************************************************
 */

int
Fire_FuelCombustion (FuelCatalogPtr catalog, size_t model )
    //FuelCatalogPtr catalog;     /* FuelCatalogData instance pointer           */
    //size_t         model;       /* fuel model id number         [0-maxModels] */
{
    size_t particle, size, life;

    double sizeClassAreaWtg[FIRE_LIFE_CATS][FIRE_SIZE_CLASSES];
    double lifeLoad[FIRE_LIFE_CATS];
    double lifeArea[FIRE_LIFE_CATS];
    double lifeSavr[FIRE_LIFE_CATS];
    double lifeHeat[FIRE_LIFE_CATS];
    double lifeSeff[FIRE_LIFE_CATS];
    double lifeEtaS[FIRE_LIFE_CATS];

    double totalArea;
    double fineLive;
    double beta;
    double betaOpt;
    double sigma;
    double ratio;
    double aa;
    double sigma15;
    double gammaMax;
    double gamma;
    double c;
    double e;

    /* Validate catalog and fuel model existence. */
    assert(catalog!= NULL && FuelCat_MagicCookie(catalog)==FIRE_CATALOG_MAGIC);
    if ( ! Fire_FuelModelExists(catalog,model) )
    {
        sprintf(FuelCat_Error(catalog),
            "Fire_FuelCombustion(): el modelo de combustible %d no existe en el cat醠ogo de combuestibles \"%s\".",
            model, FuelCat_Name(catalog));
        return (FuelCat_Status(catalog) = FIRE_STATUS_ERROR);
    }

    /* Initialize the model's fuel particle dependent variables. */
    for ( particle=0; particle<Fuel_Particles(catalog,model); particle++ )
    {
        Fuel_AreaWtg(catalog,model,particle)     = 0.;
        Fuel_SizeAreaWtg(catalog,model,particle) = 0.;
        Fuel_Moisture(catalog,model,particle)    = 0.;
    }

    /* Initialize the model's fuel combustion variables. */
    /* The following are calculated by this function. */
    Fuel_FineDead(catalog,model)        = 0.0;
    Fuel_LiveMextFactor(catalog,model)  = 0.0;
    Fuel_BulkDensity(catalog,model)     = 0.0;
    Fuel_ResidenceTime(catalog,model)   = 0.0;
    Fuel_PropFlux(catalog,model)        = 0.0;
    Fuel_SlopeK(catalog,model)          = 0.0;
    Fuel_WindB(catalog,model)           = 0.0;
    Fuel_WindE(catalog,model)           = 0.0;
    Fuel_WindK(catalog,model)           = 0.0;

    for (life=0; life<FIRE_LIFE_CATS; life++)
    {
        Fuel_LifeAreaWtg(catalog,model,life) = 0.;
        Fuel_LifeRxFactor(catalog,model,life) = 0.;
        lifeLoad[life] = 0.;
        lifeArea[life] = 0.;
        lifeSavr[life] = 0.;
        lifeHeat[life] = 0.;
        lifeEtaS[life] = 0.;
        lifeSeff[life] = 0.;
        for ( size=0; size<FIRE_SIZE_CLASSES; size++ )
            sizeClassAreaWtg[life][size] = 0.;
    }

    /* Initialize the model's fire behavior variables. */
    /* These are calculated by Fire_SpreadNoWindNoSlope(). */
    Fuel_Spread0(catalog,model)         = 0.;
    Fuel_RxIntensity(catalog,model)     = 0.;
    Fuel_HeatPerUnitArea(catalog,model) = 0.;

    /* Initialize the model's fire behavior variables. */
    /* These are calculated by Fire_SpreadWindSlopeMax(). */
    Fuel_SpreadMax(catalog,model)       = 0.;
    Fuel_AzimuthMax(catalog,model)      = 0.;
    Fuel_EffectiveWind(catalog,model)   = 0.;
    Fuel_PhiSlope(catalog,model)        = 0.;
    Fuel_PhiWind(catalog,model)         = 0.;
    Fuel_PhiEffWind(catalog,model)      = 0.;
    Fuel_LwRatio(catalog,model)         = 1.;
    Fuel_Eccentricity(catalog,model)    = 0.;
    Fuel_WindLimit(catalog,model)       = 0;

    /* Initialize the model's fire behavior variables. */
    /* These are calculated by Fire_SpreadAtAzimuth(). */
    Fuel_SpreadAny(catalog,model)       = 0.;
    Fuel_AzimuthAny(catalog,model)      = 0.;
    Fuel_ByramsIntensity(catalog,model) = 0.;
    Fuel_FlameLength(catalog,model)     = 0.;
    Fuel_ScorchHeight(catalog,model)    = 0.;

    /* Initialize the model's environmental variables. */
    Fuel_WindSpeed(catalog,model) = 0.;
    Fuel_WindDir(catalog,model)   = 0.;
    Fuel_Slope(catalog,model)     = 0.;
    Fuel_Aspect(catalog,model)    = 0.;
    for ( size=0; size<FIRE_MCLASSES; size++ )
        Fuel_EnvMoisture(catalog,model,size) = 0.;

    /* Initialize the model's combustion flag. */
    Fuel_CombustionFlag(catalog,model) = 1;

    /* If the model has no particles, we're all done. */
    if ( Fuel_Particles(catalog,model) <= 0 )
        return (FuelCat_Status(catalog) = FIRE_STATUS_OK);

    /* Initialize local fuel bed and combustion variables. */
    beta = betaOpt = sigma = ratio = aa = sigma15 = 0.;
    gamma = gammaMax = c = e = fineLive = totalArea = 0.;

    /* Accumulate surface areas by life category for the entire fuel bed. */
    for ( particle=0; particle<Fuel_Particles(catalog,model); particle++ )
    {
        life = Fuel_Live(catalog,model,particle);
        lifeArea[life] += Fuel_SurfaceArea(catalog,model,particle);
        totalArea      += Fuel_SurfaceArea(catalog,model,particle);
    }

    /* If no surface area, we're all done. */
    if ( totalArea <= Smidgen )
        return (FuelCat_Status(catalog) = FIRE_STATUS_OK);

    /* Surface area wtg factor for each particle within its life category */
    /* and within its size class category (used to weight loading). */
    for ( particle=0; particle<Fuel_Particles(catalog,model); particle++ )
    {
        life = Fuel_Live(catalog,model,particle);
        if ( lifeArea[life] > Smidgen )
        {
            Fuel_AreaWtg(catalog,model,particle) =
                Fuel_SurfaceArea(catalog,model,particle) / lifeArea[life];

            size = Fuel_SizeClass(catalog,model,particle);
            sizeClassAreaWtg[life][size] +=
                Fuel_AreaWtg(catalog,model,particle);
        }
    }

    /* Assign size class surface area weights to each particle. */
    for ( particle=0; particle<Fuel_Particles(catalog,model); particle++ )
    {
        life = Fuel_Live(catalog,model,particle);
        size = Fuel_SizeClass(catalog,model,particle);
        Fuel_SizeAreaWtg(catalog,model,particle) = sizeClassAreaWtg[life][size];
    }

    /* Derive life category surface area weighting factors. */
    for ( life=0; life<FIRE_LIFE_CATS; life++ )
        Fuel_LifeAreaWtg(catalog,model,life) = lifeArea[life] / totalArea;

    /* Accumulate life category weighted load, heat, savr, and seff. */
    for ( particle=0; particle<Fuel_Particles(catalog,model); particle++ )
    {
        life = Fuel_Live(catalog,model,particle);

        lifeLoad[life] += Fuel_SizeAreaWtg(catalog,model,particle)
                        * Fuel_Load(catalog,model,particle)
                        * (1. - Fuel_SiTotal(catalog,model,particle));

        lifeSavr[life] += Fuel_AreaWtg(catalog,model,particle)
                        * Fuel_Savr(catalog,model,particle);

        lifeHeat[life] += Fuel_AreaWtg(catalog,model,particle)
                        * Fuel_Heat(catalog,model,particle);

        lifeSeff[life] += Fuel_AreaWtg(catalog,model,particle)
                        * Fuel_SiEffective(catalog,model,particle);

        Fuel_BulkDensity(catalog,model) += Fuel_Load(catalog,model,particle);

        if ( Fuel_Density(catalog,model,particle) > Smidgen )
            beta += Fuel_Load(catalog,model,particle)
                  / Fuel_Density(catalog,model,particle);
    }

    /* Accumulate life category contribution to reaction intensity. */
    for ( life=0; life<FIRE_LIFE_CATS; life++ )
    {
        sigma += Fuel_LifeAreaWtg(catalog,model,life) * lifeSavr[life];

        lifeEtaS[life] = 1.;
        if (lifeSeff[life] > 0.)
        {
            if ( (lifeEtaS[life] = 0.174 / pow(lifeSeff[life], 0.19)) > 1.0 )
                lifeEtaS[life] = 1.0;
        }

        Fuel_LifeRxFactor(catalog,model,life) =
            lifeLoad[life] * lifeHeat[life] * lifeEtaS[life];
    }

    /* Fuel model residence time */
    Fuel_ResidenceTime(catalog,model) = 384. / sigma;

    /* Fuel model bulk density */
    if ( Fuel_Depth(catalog,model) > Smidgen )
    {
        Fuel_BulkDensity(catalog,model) /= Fuel_Depth(catalog,model);
        beta /= Fuel_Depth(catalog,model);
    }

    /* Propagating flux depends upon sigma and beta only. */
    Fuel_PropFlux(catalog,model) =
        exp((0.792 + 0.681*sqrt(sigma)) * (beta+0.1)) / (192.+0.2595*sigma);

    /* Gamma */
    betaOpt   = 3.348 / (pow(sigma, 0.8189));
    ratio     = beta / betaOpt;
    aa        = 133. / (pow(sigma, 0.7913));
    sigma15   = pow(sigma, 1.5);
    gammaMax  = sigma15 / (495. + 0.0594*sigma15);
    gamma     = gammaMax * pow(ratio, aa) * exp(aa * (1.-ratio));

    /* Factor gamma into life category reaction intensity contribution. */
    for ( life=0; life<FIRE_LIFE_CATS; life++ )
        Fuel_LifeRxFactor(catalog,model,life) *= gamma;

    /* Slope and wind intermediates constants for the fuel model. */
    Fuel_SlopeK(catalog,model) = 5.275 * pow(beta, -0.3);