• BME 200
  • Mass Balance Model of Ventilation

    System definition

    The lungs comprise a branching structure that starts off with the trachea and ends with 300 million tiny sacs called alveoli. The branching of the airways occurs over 23 generations. At generation 17 and below (called the respiratory zone) oxygen is added to the blood and carbon dioxide is removed. If we wanted to model the mass (or volume) flow of oxygen into and out of the respiratory zone, we might draw the system and CV shown below.

    Mass Balance

    Given this system, we now want to write the mass balance for oxygen in the system. We start with the general equation

    $$\frac{dm}{dt} = \sum_{i} \dot{m}_{i,in} - \sum_{i} \dot{m}_{i,out} + \sum_{i} \dot{m}_{i,gen} - \sum_{i} \dot{m}_{i,con}$$

    We assume:

    Given these assumptions, we have the mass balance

    $$0 = \dot{m}_{in} - \dot{m}_{out} + 0 - 0$$

    Since physiologists usually work with volumes instead of masses, we convert our mass balance to a volume balance because we assume that the density of oxygen gas is constant. Rewriting the mass balance we have

    $$0 = \dot{V}_{in} - \dot{V}_{out} + 0 - 0$$

    Oxygen enters the lungs during inhalation and leaves via exhalation. The volume of air (tidal volume) that enters the lungs is the same during inhalation and exhalation (note–this is the total volume of gas, not just the volume of oxygen). This is called the ventilation rate and can be denoted \(\dot{V}_A\). Since only a fraction of the inhaled and exhaled gas consists of oxygen, we multiply this fraction by the total air volume inhaled and exhaled to get the amount of oxygen. This fraction is called the partial pressure and explained in more detail here. The fraction of oxygen to air can be represented by \(\left (\frac{P_{O2}}{P_T}\right )\) where \(P_{O2}\) is the partial pressure of oxygen and \(P_T\) is total air pressure. Since there is more water vapor in the alveolar space than in the atmostphere, the partial pressure of oxygen will change slightly. We can make this distinction by denoting atmospheric and alveolar oxygen partial pressure with \(\left (\frac{P_{O2}}{P_T}\right )_{atm}\) and \(\left (\frac{P_{O2}}{P_T}\right )_{alv}\), respectively. Oxygen is also lost through diffusion through the alveoli into the blood. We denote this \(\dot{V}_{O2}\). Rewriting the volume balance to account for these details gives

    $$0 = \dot{V}_{A}\left (\frac{P_{O2}}{P_T}\right )_{atm} - \dot{V}_{A}\left (\frac{P_{O2}}{P_T}\right )_{alv} - \dot{V}_{O2}$$

    Using algebra to rearrange the equation gives

    $$(P_{O2})_{alv} = (P_{O2})_{atm} - \left (\frac{P_T}{\dot{V}_A}\right )_{alv}\dot{V}_{O2}$$

    MATLAB

    Now that we have this equation, we can program it in to MATLAB and explore what affect oxygen consuption rates have on the ventilation rate. One version of MATLAB code for the above equation is

    % use hold on to add more than one thing on a graph
    hold on;   
    
    % oxygen partial pressure in the atmosphere
    % units: mmHg
    atmPartPres = 159;
    
    % oxygen consumption rate by the body
    % units: L/min
    oxyConRate = 0.250;
    
    % total atmospheric pressure
    % units: mmHg
    totPres = 760;
    
    % alveolar ventilation rate
    % create a variable with range from 1 to 20
    % units: L/min
    alvVent = 1:20;
    
    alvParPres = atmPartPres - (totPres./alvVent)*oxyConRate;
    
    % get rid of any partial pressures less than zero
    for i = 1:20
        if alvParPres(i) < 0
            alvParPres(i) = 0;
        end
    end
    
    plot(alvVent,alvParPres,'o'); % plot(xvalue, yvalue)
    
    xlabel('Alveolar ventilation (L/min)'); % name of x-axis and units
    ylabel('Alveolar oyxgen partial pressure (mmHg)') % name of y-axis and units
    
    % change oxygen consumption rate by the body
    % units: L/min
    oxyConRate = 1.0;
    
    % calculate new partial pressures from the updated consumption rate
    alvParPres = atmPartPres - (totPres./alvVent)*oxyConRate;
    
    % get rid of any partial pressures less than zero
    for i = 1:20
        if alvParPres(i) < 0
            alvParPres(i) = 0;
        end
    end
    
    % add to plot but with different symbols
    plot(alvVent,alvParPres,'+'); % plot(xvalue, yvalue)
    

    The above code produces the following output





    Last updated:
    August 26, 2018