Estimation Problem E1-1:

Of all the N2 molecules that you take in at each breath, how many of those N2 molecules had
Galileo breathed in his entire life?

kg of (Subscript[N, 2]) in atmosphere (wikipedia), noting that 78% of atmosphere is nitrogen. This doesn’t account for the variation in concentration as a function of altitude. Note also that the 78% concentration is in terms of volume and not mass. (The mole mass of O is 16 and that of N is 14 and to be precise, we would have to correct the volume percentages for the different masses.)

kgN2Atmosphere = 5 10^18 .78

Using molecular weight of (Subscript[N, 2]) (28) this corresponds to how many moles of  (Subscript[N, 2]) molecules in atmosphere?

molesN2Atmosphere =kgN2Atmosphere/ 28

According to wikipedia, there is 0.5 liters air at STP per breath:
22.4 liters of air is one mole of ideal gas , thus
each breath is 0.5/22.4 moles, of which 0.78 is (Subscript[N, 2])
On average, one takes about 18 breaths per minute

molesN2PerBreath =0.78 0.5/22.4

How many breaths did Galileo take?
Galileo lived from 1564 to 1642

galieoLifeSpan =1642-1564

galileoBreathsLifetime = galieoLifeSpan 18 60 24 3

Assume that of all the (Subscript[N, 2]) that Galileo breathed, about 2/3 were recirculated.

molesN2GalileaLifefime = galileoBreathsLifetime molesN2PerBreath 1/3

Concentration of Galileo’s  (Subscript[N, 2]) in the entire atmosphere; this will be assumed to be completely mixed homogeneously today.

fractionGGN2Atmosphere = molesN2GalileaLifefime/molesN2Atmosphere

molesPerBreath = molesN2PerBreath fractionGGN2Atmosphere

moleculesPerBreath = 6 10^23 molesPerBreath

⟹ The estimate is 3 billion molecules per breath

Individual (Handworked) Problem I1-1: behavior of
 f (α, z) = αz(1 − z) + z log z + (1 − z) log(1 − z)

Make an assignment for expression in the assignment

theExpression = alpha x (1-x) + x  Log[x] + (1-x)Log[1-x]

It will make things simpler if we get an idea of what the function looks like.  Here I will use the Manipulate function. Change the value of alpha and see how the number of extrema changes.

Manipulate[
Plot[alpha x (1-x) + x  Log[x] + (1-x)Log[1-x],{x,0,1}],{alpha,-1,4}]

To find the extrema,  the derivatives should be computed

derivativeX= D[theExpression,x]

Let alpha = 5/2, numerically solve for the x that minimizes the solution

Plot[(theExpression/.alpha-> 5/2),{x,0,1}]

Plot[(derivativeX/.alpha-> 5/2),{x,0,1}]

From the above graph, there appears to be three local extrema, plus the two at the edges of the domain (x=0, and x=1)

The extrema at the edges of the domain will not be affected by the value of alpha, so they can be obtained immediately.

theExpression/.x-> 0
theExpression/.x-> 1

However, the Log[0]’s will make the expression indeterminant. Thus, the limits must be evaluated

Limit[theExpression,x-> 0]

Limit[theExpression,x-> 1]

Thus, for x=0 and x=1, the limiting value is 0.  Whether these are local maxima or minima depends on whether the function is increasing or decreasing at the boundary.

Find the roots to this equation numerically, starting from the left root and working to the right

We start near the point x=0.2, which is near an extremum as we have seen from the plot

tempRoot = FindRoot[(derivativeX/.alpha-> 5/2) == 0,{x,.2}]

We plug the x-value back into the function in order to obtain the y-coordinate:

(theExpression/.alpha-> 5/2)/.tempRoot

The graph shows the second local extremum near x=0.5 and it is in fact exactly at that location:

tempRoot = FindRoot[(derivativeX/.alpha-> 5/2) == 0,{x,.5}]

Corresponding y-value:

(theExpression/.alpha-> 5/2)/.tempRoot

Third one near x=0.7

tempRoot = FindRoot[(derivativeX/.alpha-> 5/2) == 0,{x,.8}]

Corresponding y-value:

(theExpression/.alpha-> 5/2)/.tempRoot

⟹, there are three maxima (x=0,x=1/2,x=1) corresponding to y=0,-0.07,0 and two minima at x=.14 and x=.86 corresponding to -.1

It seems we repeated the above steps over and over again, why not generalize as a function

aRoot[anAlpha_,near_] := FindRoot[(derivativeX/.alpha-> anAlpha) == 0,{x,near}]  (* find the root close to the x-value 'near' *)

theValueAtRoot[anAlpha_, near_] := (theExpression/.alpha-> anAlpha)/.aRoot[anAlpha,near]  (* get the corresponding y-value *)

Employing the functions for α=3/2

Plot[(theExpression/.alpha-> 3/2),{x,0,1}]

Plot[(derivativeX/.alpha-> 3/2),{x,0,1}]

aRoot[3/2,.5]

theValueAtRoot[3/2,.5]

⟹, there are two maxima (0;0) and(1;0) and one minimum at (0.5;-0.32)

The function’s character changes when the two inflection points converge to 1/2. Thus, we can use the second derivative to find the α where this will happen

secondDerivativeX= D[theExpression,{x,2}]

Here are the solutions for the positions of the inflection points:

xSols =Solve[secondDerivativeX==0,x]

The two solutions converge when the discriminant vanishes. Thus the critical value of alpha is

Solve[-2 alpha+(alpha^2)==0, alpha]

⟹ α = 2 is where the number of maxima changes.

For plotting the positions of the inflection points, a parametric plot is useful

ParametricPlot[
{
{x/.xSols[[1]],alpha},
{x/.xSols[[2]],alpha}
}
,{alpha,0.1,2.9},AxesLabel-> {"x","α"},PlotLabel-> "α(x=inflection point)", AspectRatio-> 1]

Individual Problem I1-2: 

Factor x^5 − 11x^4 + 41x^3 − 61x^2 + 30x

Factor[x^5 - 11 x^4 + 41 x^3 - 61 x^2 + 30 x]

This fifth order polynomial has zeros at x=0,1,2,3,and 5

Individual Problem I1-3: 

∫(e^(x^x))(x^x)(log(x)+1)cos((e^(x^x))+1) dx

define expressions for the integrand and integral

integrand = Exp[x^x] x^x (Log[x] + 1) Cos[Exp[x^x] + 1]

integral = Integrate[integrand,x]

It is not obvious from the above what Mathematica is doing with the constant of integration, we could figure this out from a known case; these will show that Mathematica sets the integration constant to zero.

Integrate[x^2,x]
Integrate[x^2,{x,0,y}]
Integrate[Sin[x],x]

⟹, Mathematica sets the integration constant to zero.

Plot[{integrand,integral}, {x,0,1},PlotStyle-> {{Thick, Red},{Thick, Blue}},
PlotLabel->Row[ {Style["Integrand",{Red}],Style["Integral",{Blue}]},"\t"]]

The above shows quite nicely that the integral’s maximum corresponds to the integrand’s zero, which is related to the first fundamental theorem of calculus.

Finding the roots to the integral and the integrand

Again we have to start with some initial values (0.2, 0.35 and 0.7). We read the off the graph above.

FindRoot[integrand==0,{x,0.35}]

FindRoot[integral==0,{x,0.2}]

FindRoot[integral==0,{x,0.7}]

⟹ The integral has roots at 0.137 and 0.663, and then integrand at 0.377

Individual Problem I1-4: 

First one hundred primes, and histogram of their differences

The list of one hunderd primes

oneHundredPrimes = Table[Prime[i],{i,1,100}]

A plot of the first 1000 primes.

ListPlot[Table[Prime[i],{i,1,1000}],AxesLabel-> {"n","\!\(\*SuperscriptBox[\(n\), \(\(th\)\(\\ \)\)]\)Prime"},ImageSize-> Large, BaseStyle-> Large]

Histogram of the difference of the first 1000 primes.

Histogram[Differences[Table[Prime[i],{i,1,1000}]],25,FrameLabel-> {"Difference betweeen Primes", "Number"},Frame-> True, PlotLabel-> "Histogram of differences of the first 1000 primes",BaseStyle-> Large, ImageSize-> Full]

Challenge, counting instances of differences not equal to 2

A function to construct the differences list

primeDifferences[n_] := Differences[Table[Prime[i],{i,1,n}]]

example

primeDifferences[12]

A function to select only those differences that are prime

primeDifferenceAlsoPrime[n_] := Cases[primeDifferences[n],_?PrimeQ]

example

primeDifferenceAlsoPrime[12]

countPrimes[n_] := Count[primeDifferenceAlsoPrime[n],a_/;a≠2]

example

countPrimes[1000]

⟹ The only difference which occurs between the prime numbers which is again a prime number is the number 2.

Individual  Problem I1-5: 

code snippets:
i = 100; Table[i++, {i, 1, 20}]
j = 1; Table[j++, {i, 1, 20}]

(a) fOne[i] = {i, DateString[]} Table[Pause[1];fOne[j],{j,1,5}]
(b) fTwo[i_] = {i, DateString[]} Table[Pause[1]; fTwo[j], {j, 1, 5}]
(c) fThree[i_] := {i, DateString[]} Table[Pause[1]; fThree[j], {j, 1, 5}]

part 1

i = 100; Table[i++, {i, 1, 20}]

i = 100; Table[i++, {i, 1, 20}]

Here, the initial assigment to i has no effect on the output, the second command just prints a table and uses the internal values of i. The ++ also has no effect, since it is a post-increment. The command i++ FIRST returns i and THEN increments it. The command ++i FIRST increments the number and THEN returns the incremented value.

j = 1; Table[j++, {i, 1, 20}]

j = 1; Table[j++, {i, 1, 20}]

Here,  the value of j is being updated, and the Table’s internal variable i is just counting the steps.  To see what the difference is, one could do the following experiment:

j = 100; Table[j++, {i, 1, 20}]

part 2

fOne[i] = {i, DateString[]} 
Table[Pause[1];fOne[j],{j,1,5}]

fOne[i] = {i, DateString[]} 
Table[Pause[1];fOne[j],{j,1,5}]

Here, the first line only initializes f[100] (since the previous of i was 100) and then displays the result. The next line accesses fOne[1] to fOne[5], but the system had not initialized their valie; so mathematica just returned all it knew about the function.

(b) fTwo[i_] = {i, DateString[]} 
Table[Pause[1]; fTwo[j], {j, 1, 5}]

fTwo[i_] = {i, DateString[]} 
Table[Pause[1]; fTwo[j], {j, 1, 5}]

Here, fTwo was immediately defined for any pattern i_ and the right hand side evaluated to i (which was 100) and the current date at the time of the evaluation.  So, f[ anything] would always have the same value in the future.

fThree[i_] := {i, DateString[]} Table[Pause[1]; fThree[j], {j, 1, 5}]

fThree[i_] := {i, DateString[]} 
Table[Pause[1]; fThree[j], {j, 1, 5}]

Here:
1) The result of the first line was never displayed, because nothing was set or computed when the evaluation took place.  
2) A “set delayed “:=” was used which means, “define the function, but apply its definition when I use it in the future.
3) Thus, each member of the list is a “j” according the value of “j” when the function was evaluated, and the time at which it was evaluated.

Individual Problem I1-6: 

Eigenvalues and eigenvectors of  (({{7, 6}, {6, 2}}))

theMatrix = {
{7,6},
{6,2}
}

MatrixForm[theMatrix]

Eigenvalues[theMatrix]

The eigenvectors are given in list form:

eVects = Eigenvectors[theMatrix]

MatrixForm[eVects[[1]]]

MatrixForm[eVects[[2]]]

We may not know what eigenvalues and eigenvectors are right now, but we know we can compute them if we need them in the future.

Individual Problem I1-7: 

solutions to the ordinary differential equation:
({{ ((dy(x))/(dx)) ,  +y(x)==sin(x) }, { y(x==0)==b ,  }})

Set up assigments for use in DSolve

lhsODE = D[y[x],x] + y[x]

rhsODE = Sin[x]

Set the solutions equal to a variable called dsols. There is only one solutions because our equation is first order and linear.

dsols = DSolve[{lhsODE==rhsODE,y[0]==b},y[x],x]

theSol = y[x]/.dsols[[1]]

Manipulate[
Plot[theSol/.b-> beta,{x,0,2}],{beta, -1, 1}]

Plot[
Evaluate[
Table[theSol/.b-> beta,{beta, -1, 1,.2}]],{x,0,7}]