[boost] Re: physical quantities (was: runtime unit-modifiable ...)

From: Andy Little (andy_at_[hidden])
Date: 2004-01-07 12:28:29

"David B. Held" <dheld_at_[hidden]> wrote
> "Andy Little" <andy_at_[hidden]> wrote

> > [...]
> > Feel free to recode  in your 'problem domain'  with m and s my
> > friend. I'll take my syntax .
>
> But you didn't show the code using your syntax, so there's nothing
> to compare.  You show the code with your syntax, and I'll show it
> with postfix units, and then we'll have apples and apples.

code  redone with pqs-100-01 syntax below  (compile checked only)

regards
Andy Little

/*
quick demo of previous post
redone with physical-quantity type
note that angle are incomplete
compiles on VC7.1 not tested on gcc3.2
not run checked(will give div 0) */

#include "../physical_quantities/pqs.hpp"
using namespace physical_quantities;

// class for info abour rotor
struct RotorDialog{
    q_length::mm            m_outer_dia;
    q_velocity::m_div_s     m_vw;
    double                  m_tsr;
    double                  m_ellipticality;
    q_density::kg_div_m3    m_rho_air;
    double                  m_cd;
    double                  m_cl;
    int                     m_numblades;
    //should be q_angle
    double                  m_alpha;
}rotordialog;

// lash-up to return chord and angle
struct chord_omega{
// 2nd param should be q_angle::rads
    chord_omega(q_length::m,double){}
};

// some func
chord_omega getChordOmega(q_length::m const& r);

const double PI = 3.1414;

// section spacing
inline q_length::mm  get_dr()
{
    return q_length::mm(100);
}

int main()
{
    getChordOmega(q_length::m(1));
}

//redo of the original function
//some now unnecessary scalings removed

chord_omega getChordOmega(q_length::m const& r)
{
    const q_length::m tip_r
    = rotordialog.m_outer_dia /2;
    const q_length::m rr
    = r;
    const double& tsr
    = rotordialog.m_tsr;
    const double& ellipticality
    = rotordialog.m_ellipticality;
    const q_velocity::m_div_s Vin_y
    = rotordialog.m_vw;
    const q_length::m  dr
    = get_dr();
    const q_time::s dt(1);
    const q_velocity::m_div_s Vout_y
    = Vin_y
    * ( (rr < tip_r)
        ? 1.0L
            - ((2.0L / 3.0L)
                * (1.0L - ellipticality)
                * (ellipticality
                    * sqrt(1.0L - ( power<2>(rr) )
                                    /( power<2>(tip_r) ))))
        : 1.0L - ((2.0L/3.0L)
            * (1.0L - ellipticality)) );
    const q_velocity::m_div_s Vb_y
    = (Vin_y + Vout_y) / 2.0L;
    //gives Dimensional Analysis error///
    //const q_mass::kg mass
    //= Vb_y * 2 * PI * (rr)
    //* rotordialog.m_rho_air * dr;
    /////////////////////////////////////
    const q_mass::kg mass
    = Vb_y * 2 * PI * rr
    * rotordialog.m_rho_air * dr * dt;
    const q_force::N Fb_y
    = mass * (Vin_y - Vout_y) / dt;

    q_velocity::m_div_s  Vout_x(0);
    q_velocity::m_div_s oldVout_x(0);
    const double cl = rotordialog.m_cl;
    const double drag_coeff
    = rotordialog.m_cd / cl;
    for (int i = 0 ;i < 1000 ; ++i){
        const q_velocity::m_div_s Vb_x
        = Vout_x / 2.0L
        + Vin_y * tsr * (rr / tip_r);
        // note beta should be q_angle::rads
        const double beta
        = atan2(Vb_y(),Vb_x());
        const double sinbeta
        = sin(beta);
        const double cosbeta = cos(beta);
        const q_force::N Lift
        = Fb_y / (cosbeta + drag_coeff * sinbeta);
        const q_force::N Fb_x
        = Lift * (sinbeta - drag_coeff * cosbeta);
        oldVout_x = Vout_x;
        ///////// gives  D.A error
        //Vout_x = Fb_x / mass;
        ///////////////////////
        Vout_x = (Fb_x / mass) * dt;
        /*abs(Vout_x - oldVout_x)...no abs()*/
        if ( (Vout_x - oldVout_x) < 0.0001L * Vin_y ){
            q_length::m chord
            = Lift
            / ((power<2>(Vb_x) + power<2>(Vb_y) )
                * 0.5 * rotordialog.m_rho_air
                * cl * dr * rotordialog.m_numblades);
            return chord_omega(chord,
                        beta - rotordialog.m_alpha
                        * 0.017453293L);
        }
    }
    return chord_omega(q_length::mm(0),0);//fail
}