INTRODUCTION TO MATHCAD
Mathcad is a "numerical scratchpad" which allows you to solve complex systems of equations, find roots of polynomials and plot results.  It can evaluate integrals and derivatives, both numerically and symbolically.  Units are easily incorporated, and numerical answers cannot be found if the units are inconsistent.  Mathcad is a powerful tool, but easy to use.  These notes are a brief introduction.  They are meant to be used while working in Mathcad.  

If you enter Mathcad, you will notice there are several menus that appear on the top line (File, Edit, Text, etc.).  There is a drop-down "Math" palette containing icons for calculator, graphs, matrices, etc.   If you move the cursor to an icon, there is a brief message explaining it.  If you click on an icon, a palette with various operations appears.  Click on the calculator palette and select the  symbol.  Type = and you get 

Some symbols, such as , g (acceleration of gravity) and R (degrees Rankine), have a specific meaning unless you redefine them.
Now, choose the nth root symbol and put in some numbers, such as:

To enter a number, just click on the appropriate placeholder (the black box).

These examples show how to use the "=" sign, which means evaluate or calculate.  Mathcad can be used as a calculator in this way, although it is a bit expensive and not particularly portable (unless you have a lap top).  If you choose the Evaluation and Boolean icon, you will see there are several kinds of equals, in addition to inequalities.  The ":=" sign means "is defined by."  Type x followed by this sign, and assign a number such as:

In addition to using the Evaluation palette, this ":=" sign can be inserted simply by typing a colon (:).  Once x is defined, we could then use it in an expression such as:

To get the natural log, simply type "ln".   There are a number of built in functions available as you can see by clicking of the "f(x)" icon and scrolling through the list.  It includes error functions and Bessel functions, in addition to the more mundane ones.  Each comes with a brief explanation.  There are also functions for numerical integration using the Runge-Kutta method.  We will use the simple ones here.

Return to the present problem.  If you now type a regular equals sign, from the palette or the keyboard, you get:

You could also define a dependent variable, f, as this function of x, for use in later calculations:

If you want the numerical value of f, simply type 

This type of function may be called explicit:  you can write f equals some expression which contains only the independent variable, x.  Suppose, however, we turned the equation around and wanted to find the x for which:

This is an example of an implicit equation.  You cannot rearrange this to get it into the form x = (something containing no x).  You can, however, solve this equation in Mathcad.  This is accomplished
 in a "given - find" block.  Mathcad does this numerically and needs initial guesses for the unknowns, in this case x.    Enter a guess:

Then type "given" with or without a capital G.  (In most cases, Mathcad is not case sensitive.)

Next, enter the equation.  In these blocks you must use a third type of "equals" which means "is constrained by" and is used in these "given-find" blocks.  (In the Evaluation palette, this is called a "Boolean equals.")  This is the bold faced "=" from the evaluation palette (or CTRL = from the keyboard).  

To get an answer, type "find(x)" followed by an =

If you next type x=, you will not see the solution (4.75), but the initial guess.  This is because you have not redefined the value of x.  If you want the value of x to be the solution, you can proceed as follows:

(Initial guess)

Now x has been replaced by the solution of the equation.  To see the value, just type:

The major restriction on the given-find statements is that you cannot use general subscripted variables (i.e., i or j) in them.  You can, however, solve some interesting problems.  In the analysis of flow of a power law fluid in an annulus, you need to find the value of  in the following equation:

where  is the ratio of the inner to the outer radius, n is the power law index and  is the radial coordinate.   In order to get the integral sign, go to the calculus palette and choose the definite integral.  We can solve this equation using "given - find" as above.  First define the  and n and guess a value for :

We can also include units in Mathcad calculations.  Many units are already defined in Mathcad;  others can be assigned.  If you click on the measuring cup, you can see a list of the units available.  You can change the default system of units by choosing Math then Units then Change System of Units.  However, whatever system of units you are in, you can also always go from one to another.  For example, take the acceleration of gravity, g.  Type g= to see its value in SI (presuming the default system is set on MKS):

You notice a "placeholder" or a black square after the units.  Click on this value and insert another set of units:

As an example of a problem with units, suppose we want investigate the P-V-T behavior  of methane using the Redlich-Kwong equation of state.  First, we enter the critical values for methane:

When you do this, the pressure units, bar, are highlighted.  Mathcad doesn't recognize that unit.  We can define this as a new unit using the third type of equal sign:

Here is yet another type of equals.  You can get this symbol from the Evaluation Palette, or by typing a tilde (~).  Once you define the units bar in this fashion, the highlight disappears.   Notice, however, what happens if you use the ":=" to define the bar unit.  In this case, the bar in the Pc equation still remains highlighted.  This "triple equals"  is a "global assignment;"  it's good throughout the Mathcad document.  The other equals (:=) means the definition is good from this point on, and can be re-defined.

Now, introduce the gas constant and the two constants associated with the Redlich-Kwong equation:

We can define the moles in terms of grams for methane:
(See what happens later if you do not use this global 
assignment.) 

Suppose we were given a volume and a temperature and wanted to find the pressure.  The specified temperature and volume are:

According to the Redlich - Kwong equation, pressure is given by the expression:

Therefore:

To find the volume at 165K and 19.39 bar, we have to solve an implicit equation using Given - Find: 

(Initial guess)

The default units of kg and m can be changed to moles and cm, or anything else you desire.  

In this problem, as in many others, the solution is sensitive to the initial guess.  Try using a different value for the guess.  Sometimes, the numerical procedure will not converge, or will converge on a different root. 

One of the other Mathcad features is root finding.  Given a polynomial, we can find all the roots, real and imaginary, using the "polyroots" function.  To use this feature in this problem, first, write the equation in polynomial form.  The Redlich - Kwong equation is cubic in the volume (it is an example of a "cubic equation of state"):

Next, form a vector, call it V, of the coefficients in this equation. The first element is the coefficient of V0, the next, V1, etc:

Notice that the elements of this vector, as with any vector, must have the same units.  To use the root finding function, however, the elements must all be dimensionless.  The three roots are then found simply by using the "polyroots" function:

To decide on which root means what, we have to know a little about the physics of the situation.
The vapor volume is the maximum value, so define:

Since we took out the units in the vector V, we have 
to add them back again.  We also have to put back 
the same ones we took out, that is, if we multiplied the
elements of V by kg/m3, we must use m3/kga s the units for V.

The liquid volume is the minimum:

The middle value has no physical significance.  (It might be interesting to some, at least, to note that the actual vapor volume of methane at these conditions is 0.032, as calculated.  The liquid volume is 3.095 10-3.)

Try using a different T and P, say 300K and 1 atm, and see what the roots are.
Now, let's plot the compressibility factor for methane as a function of pressure for a temperature of 300K.  We need to find v for a number of values of P.  However, only one given-find statement is needed.  We can proceed as before:

(We already have an initial guess for v above.)

The v(P) means that v is a function of P;  for each value of P, a different v is found. If we specify a particular value of P we get a numerical answer:

Now, introduce a range of pressures.  This could have been done before the given - find statement if desired:

This means P takes on the values 0.1, 0.2, 0.3, ... all the way to 10 MPa.  The two dots appear when you type a semicolon (;).  Note that units must be included.
Define the compressibility factor:

Notice that z is also a function of P.  Now create an x-y plot using the Graphics menu,  the Graphing palette, or simply SHFT @.  Put the cursor on the placeholders for the x and y axes and put in P and z(P), respectively:

Notice the values of P are much larger than you defined. That's because this P is in Pa (the default units for pressure), not MPa.  To change the scale, you can plot P/MPa rather than P (or P/atm, P/psi, or anything you want).  If you double click on the graph, or if you choose "Format" with the right mouse buton, you can customize your graphs.  By clicking the mouse in the region near the graph (or anything else), you can get a line enclosing it.  You can move the graph, text or equation anywhere on the page using the cursor.  When the cursor is at the boundary of the line, an arrow appears.  This can be used to re-size the graph.
As an example of a somewhat more complicated problem, consider the velocity of a freely falling drop of water in air.  
The diameter of the drop is:

Density of the water:

Density and kinematic viscosity of air:

The drop velocity v depends on the friction factor f, which also depends on the Reynolds number, NRe, vD/.  We can write this as 2 equations with 2 unknowns, v and f.  In addition, if we are above the creeping flow regime, f is given by one of two expressions, depending on the value of the Reynolds number.  If NRe is greater than 500, f is 0.44.  If NRe is less than 500, f depends on NRe.  The "if" statement in Mathcad provides an easy way to handle this.  The syntax of this statement is:  if (condition, true, false).  In this example, f equals 0.44 if vD/ is greater than 500.
Initial guesses:

Momentum balance:

Friction factor:
The solution to this given - find statement will be a two component vector, v and f.  To form this vector, go to the Vectors and Matrices Palette, choose the symbol for a matrix, enter the number of rows and columns (2 and 1 in this case), and hit "Create".  You will then see a general 2 component vector.  Define the elements as the unknowns and define this vector as the "find":

The value of the Reynolds number is:

There is often an easier way to set a problem up, and so make the solution easier.  Very often, with  multiple unknowns in a given - solve block, the solution is very sensitive to the initial guesses.  Therefore, it is always easier to have only one unknown.  In the present problem, we are interested in the velocity, and not really the Reynolds number and friction factor.  We could, therefore, define these two functions outside of the given - find statement.  

To avoid any possible confusion between these 2 examples, we'll use a different symbol for velocity. 
Initial guess for vel:

If we wanted the Reynolds number and the friction factor, then we
could evaluate them:

While it is convenient to have only one unknown in a given - find block, it is not always possible to do so.  However, the fewer the unknowns, the easier it is to find a solution.  The number can usually be reduced to a number smaller than you might think.  It just requires a little "up-front" thought about the nature of the problem and what are the knowns.  Consider this complex vapor - liquid equilibrium problem.  A solution of 2-propanol and water is at a pressure of 101.33 bar and the mole fraction of propanol in the vapor phase is 0.4.  We want to find the temperature and the mole fraction of propanol in the liquid.   The liquid phase is non-ideal,  and we have to use activity coefficients for the liquid.  These depend on temperature and liquid composition.  We'll assume the Wilson model can be employed.  

First, introduce some constants for this system of components:

Two temperature-dependent functions which appear in the activity coefficient model are:

We also have the Antoine equations for the vapor pressures:

The mole fraction of component 2 is:

The activity coefficients are:

The total pressure is:

Notice that this differs from the case of Raoult's law becuase of the presence of the activity coefficients, 1 and   
The vapor phase composition is given by:

There are a large number of unknowns in this problem:  the vapor pressures, the activity coefficients, the liquid phase compositions, the temperature.  However, there are only two independent variables, x1 and T.  Once these are known, all the others can be found.  The given - find statement is just:
Initial guess:

If desired, we could find all the other dependent variables:

Mathcad cannot, apparently, handle functions which have general subscripts.  If, for example, you want to evaluate the vapor pressure of a set of components as a function of temperature (as in the previous problem), you might first introduce the Antoine constants.  To place a subscript on a variable, use the left bracket [ .  After you enter the first number, type a comma.  Mathcad will then start to form a vector.

You can either enter the constants times K, or redefine the values as done here.

You would like to define the vapor pressure with an equation like:

You can define variables with subscripts as we did here by typing "P." (P with a period).  This is sometimes convenient.  There is a fundamental difference between the subscript formed in this way and that formed with the left bracket.  The former ("sat" above) is simply part of the name of a variable.  The latter ("i") makes the variable a vector with components i.  In addition, the variable i must be defined separately.

Returning to the problem at hand, note that if you define the vapor pressure  in this way (using the := sign to define the function) Mathcad tells you this is not a proper name.  One way around this difficulty is to employ the "vectorize" option:  treat Psat as a vector make the right side of the equation a vector.  Notice that we did this sort of thing above when we said C:=C1*K.  Instead of redefining each of the components of C, we redefined the whole vector.

Write the Antoine equation without subscripts, but place an overbar (from the Vectors and Matrices Palette) over the entire right side:

Without the vectorize operation, Mathcad will want to divide one vector by another:  an illegal array operation. 

There is one problem with this approach.  Choose a temperature and evaluate the vapor pressures (in this equation, temperature must be in Celsius, but Mathcad, unfortunately, cannot go from C to K.  If you want to employ Celsius degrees, you must also use the symbol K and be careful):

How did we wind up with 3 pressures when we only defined 2 components?  Mathcad always assumes that the first element of a vector is numbered 0, the next 1, and so forth.  Since we did not define the constants A0, etc.,  Mathcad took them to be zero.  Hence the first vapor pressure is exp(0) or 1.  We can get around this problem by redefining the first element using the ORIGIN command.  (Case is important here.)

Now calculate the pressures:

The first element of the vector is now numbered 1, and we have only 2 components. 

BE CAREFUL WHEN YOU USE THIS "ORIGIN" STATEMENT!
As a final class of examples, consider curve fitting.  The first example is straightforward.  The following set of data gives dimensionless temperature  as a function of dimensionless distance  in a fiber spinning experiment.  Each set forms a vector.  One way to represent the vector is to use subscripts to denote each element of its elements.  Since, by default, the numbering of elements in a Mathcad vector begins with 0, we will return to that default setting:  
Here we used the regular definition :=, rather than the global, so that this numbering system takes effect from here on.

Mathcad assumes the increment is one unless you tell it otherwise.
Next, type in the data:

We expect the temperature to follow the equation:

where Nu is the Nusselt number, which we must find.  The constant, 0.1625, comes about from the way  and Nu are defined.

We can take the natural log of this equation to put it in linear form:

where  is the slope (-0.1625 Nu).  Finding the slope of a straight line is easy in Mathcad.  First, we need to define a new dependent variable as the log of :

The slope of a straight line through these points is found by using the "slope(x,y)" function, where x is the vector of independent variables :

Therefore, the Nu is:

To verify, we can plot the data (, ) and the prediction, p:

After you create the x-y plot, you can customize it to display the experimental data as points, rather than a line.  Double click on the graph, go to Format and Traces.   For the first trace (in this case the measured temperatures) choose a symbol, and "points" for Type.  To plot 2 variables on the same axis, put a box around the first one, here i, then type a comma. A new place holder will appear on the axis.  Unless you tell it otherwise, Mathcad will plot each of the dependent variables against the same independent variable.
In this case, the fit is not all that good, perhaps because the straight line fit does not have an intercept of zero.  We can check that with the "intercept" function:

Mathcad 2000 contains a series of curve fitting routines which can be found with the function key.  One is simply "line."   If we employ the line function here we get:

The intercept and the slope are given directly.
The predicted best fit using a straight line model with non-unity intercept is:

The fit is improved with , but now the initial temperature is not 1.  You might also try changing the y-axis to a log scale using X-Y Plot, Format, X-Y Axes.

There is one more fit we can try, this time forcing the intercept to be unity.
Define a new function, SSE(Nu), as the sum of the square error between the measured and predicted temperature:

We want to find the value of Nu which makes this function a minimum.  Notice it will not be zero unless the fit were exact.  To find this minimum, we can use a procedure related to the Given - Find method: 

Notice we already have a guess for Nu defined above.

Instead of "Find," we use "Minerr".  The value of Nu which gives the minimum error in the equation is:

This is quite different from the previous value.  We can define another predicted temperature:

And plot this function to compare with the data and the other two predictions:

Mathcad 2000 also provides an exponential fit function.  The form is a*(e^b*x)+c.  There are a few problems with this, not the least of which is that we do not want a fit of this form.  (It has the extra constant, c.  There does not appear to be any way to use this function without the constant.  Also, the results seem very sensitive to the initial guesses.  We will not pursue this method.
For those who prefer a more rigorous way of finding the best value of Nu, we can find the value which forces the derivative of the SSE equation with respect to the Nusselt number to be to zero:

This is quite close to what we just found using Minerr.  In general, the procedure using the derivative seems to give better results (in the sense of a smaller value of the sum of the square errors) than Minerr procedure.  However, the derivative method is more sensitive to the initial guess.
There are two aditional Mathcad functions which can be used for data fitting, linfit and genfit.   

Suppose we have the following set of x-y data.   We want to fit these data to an equation of the form:

We could use either of the methods discussed above.  In addition, we could use the linfit function.  This is useful whenever we can write the dependent variable as series of any functions, not just powers as in this example.

First, define the variables:

Then define a vector, called F here, which gives the various functions which appear in the expression for the dependentvariable.  In the present case:

The linfit function returns a vector, called S here, which contains the coefficients used to create a linear combination of the functions in F which give the best fit to the data in x and y:

The best fit to the data is then easily written in matrix notation:

Another way of representing this series is:

You can easily verify the two results are the same.
Very often, we cannot write the dependent variable as a linear combination of functions.  Suppose we wanted to fit the above y values to an equation of the form:

and want to find the best values of u0 and u1.  In this case, we can use the genfit function.  We first define a vector G(x,u).  The first element is the desired function;  the second, the derivative of the function with respect to the first constant;  and the third, the derivative with respect to the second:

Then we give a vector of guesses for the constants:

The genfit function returns the best values of the constants:

The best fit is then:

Clearly, the second function does not give a very good fit in this case.  The polynomial is clearly superior.