Solving Equations

6.7 A chain of resistors

Consider a long chain of resistors wired up like this:

reschain

All the resistors have the same resistance R. The power rail at the top
is at voltage V_+=5V. The problem is to find the voltages V1…VN at the internal points in the circuit.

3V_1 - V_2 - V_3 = V_+,
-V_1 + 4V_2 - V_3 - V_4 = V_+,
\vdots
-V_{i-2} - V_{i-1} + 4V_i - V_{i+1} - V_{i+2} = 0,
\vdots
-V_{N-3} - V_{N-2} + 4V_{N-1} - V_N = 0,
-V_{N-2} - V_{N-1} + 3V_N = 0.

Express these equations in vector form Av = w and find the values of the matrix A and the vector w. Write a program to solve for the values of the Vi when there are N=6 internal junctions with unknown voltages. (Hint: All the values of Vi should lie between zero and 5V. If they don’t, something is wrong.) Now repeat your calculation for the case where there are N=10,000 internal junctions. This part is not possible using standard tools like the solve function. You need to make use of the fact that the matrix A is banded and use the banded function from the file banded.py, discussed in Appendix E.

6.9 Glycolysis

The biochemical process of glycolysis, the breakdown of glucose in the body to release energy, can be modeled by the equations

{d x\over d t} = -x + ay + x^2 y,
{d y\over d t} = b - ay - x^2 y.

Here x and y represent concentrations of two chemicals, ADP and F6P, and a and b are positive constants. One of the important features of nonlinear linear equations like these is their stationary points, meaning values of x and y at which the derivatives of both variables become zero simultaneously, so that the variables stop changing and become constant in time. Setting the derivatives to zero above, the stationary points of our glycolysis equations are solutions of

-x + ay + x^2 y = 0,       b - ay - x^2 y = 0.

Demonstrate analytically that the solution of these equations is

x=b, y = {b\over{a+b^2}}.

Show that the equations can be rearranged to read

x = y(a+x^2), y = {b\over{a+x^2}}

and write a program to solve these for the stationary point using the relaxation method with a=1 and b=2. You should find that the method fails to converge to a solution in this case.  Find a different way to rearrange the equations such that when you apply the relaxation method again it now converges to a fixed point and gives a solution. Verify that the solution you get agrees with part a.

6.16 The Lagrange point

There is a magical point between the Earth and the Moon, called the L_1 Lagrange point, at which a satellite will orbit the Earth in perfect synchrony with the Moon, staying always in between the two. This works because the inward pull of the Earth and the outward pull of the Moon combine to create exactly the needed centripetal force that keeps the
satellite in its orbit. Here’s the setup:

l1

Assuming circular orbits, and assuming that the Earth is much more massive than either the Moon or the satellite, show that the distance r from the center of the Earth to the L_1 point satisfies

{GM\over r^2} - {Gm\over(R-r)^2} = \omega^2 r ,

where M and m are the Earth and Moon masses, G is Newton’s gravitational constant, and \omega is the angular velocity of both the Moon and the satellite.  The equation above is a fifth-order polynomial equation in r (also called a quintic equation). Such equations cannot be solved exactly in closed form, but it’s straightforward to solve them numerically. Write a
program that uses either Newton’s method or the secant method to solve for the distance r from the Earth to the L_1 point. Compute a solution accurate to at least four significant figures.

The values of the various parameters are:

G = 6.674\times10^{-11}\,\mathrm{m}^3\mathrm{kg}^{-1}\mathrm{s}^{-2}
M = 5.974\times10^{24}\,\mathrm{kg},
m = 7.348\times10^{22}\,\mathrm{kg},
R = 3.844\times10^8\,\mathrm{m},
\omega = 2.662\times10^{-6}\,\mathrm{s}^{-1}.

You will also need to choose a suitable starting value for r , or two starting values if you use the secant method.

6.18 The temperature of a light bulb

An incandescent light bulb is a simple device – it contains a filament, usually made of tungsten, heated by the flow of electricity until it becomes hot enough to radiate thermally. Essentially all of the power consumed by such a bulb is radiated as electromagnetic energy, but some of the radiation is not in the visible wavelengths, which means it is useless for lighting purposes.

Let us define the efficiency of a light bulb to be the fraction of the radiated energy that falls in the visible band. It’s a good approximation to assume that the radiation from a filament at temperature T obeys the Planck radiation law, meaning that the power radiated per unit
wavelength \lambda obeys

I(\lambda) = 2\pi A hc^2 {\lambda^{-5}\over{e^{hc \lambda k_B T}-1}} ,

where A is the surface area of the filament, h is Planck’s constant, c is the speed of light, and k_B is Boltzmann’s constant. The visible wavelengths run from \lambda_1=390\,nm to \lambda_2=750\,nm, so the total energy radiated in the visible window is\int_{\lambda_1}^{\lambda_2} I(\lambda)\>d\lambda and the total energy at all wavelengths is \int_0^\infty I(\lambda)\>d\lambda. Dividing one
expression by the other and substituting for I(\lambda) from above, weget an expression for the efficiency \eta of the light bulb thus:

$latex \eta = {\int_{\lambda_1}^{\lambda_2} \lambda^{-5}/(\e^{hc/\lambda
k_BT}-1) \>\dd\lambda\over
\int_0^\infty \lambda^{-5}/(\e^{hc/\lambda k_BT}-1) \>\dd\lambda}$,
where the leading constants and the area A have canceled out. Making the substitution x=hc/\lambda k_BT, this can also be written as

$latex \eta = {\int_{hc/\lambda_2 k_BT}^{hc/\lambda_1 k_BT} x^3/(\e^x-1) \>\dd x
\over \int_0^\infty x^3/(\e^x-1) \>\dd x}
= {15\over\pi^4}
\int_{hc/\lambda_2 k_BT}^{hc/\lambda_1 k_BT}
{x^3\over\e^x-1} \>\dd x$,
where we have made use of the known exact value of the integral in the denominator.
Write a Python function that takes a temperature T as its argument and calculates the value of \eta for that temperature from the formula above. The integral in the formula cannot be done analytically, but you can do it numerically using any method of your choice. (For instance, Gaussian quadrature with 100 sample points works fine.) Use your function to make a graph of \eta as a function of temperature between 300K and 10,000\K. You should see that there is an intermediate temperature where the efficiency is a maximum.
\Calculate the temperature of maximum efficiency of the light bulb to within 1K using golden ratio search. (Hint: An accuracy of 1K is the equivalent of a few parts in ten thousand in this case. To get this kind of accuracy in your calculation you’ll need to use values for the
fundamental constants that are suitably accurate, i.e.,~you will need
values accurate to several significant figures.)  Is it practical to run a tungsten-filament light bulb at the temperature you found? If not, why not?