Ian Hickson -- Tutor: Dr Laughton
for and while loops The following straight-forward for loop:
statement(1);
for (statement(2); expression(x--); statement(4)) {
/* if break 1 goto end of loop */
statement(3);
/* if break 2 goto end of loop */
}
/* end of loop */
statement(5);
...is equivalent to this while loop:
statement(1);
statement(2);
while (expression(x--)) {
/* if break 1 goto end of loop */
statement(3);
/* if break 2 goto end of loop */
statement(4);
}
/* end of loop */
statement(5);
Similarly, the following straight-forward do..while loop:
statement(1);
statement(2);
do {
/* if break 1 goto end of loop */
statement(3);
/* if break 2 goto end of loop */
statement(4);
/* if break 3 goto end of loop */
} while (expression(x--));
/* end of loop */
statement(5);
...is equivalent to this for loop:
statement(1);
statement(2);
/* if break 1 goto end of loop (after loop) */
statement(3);
/* if break 2 goto end of loop (after loop) */
statement(4);
/* if break 3 goto end of loop (after loop) */
for (; expression(x--); ) {
/* if break 1 goto end of loop */
statement(3);
/* if break 2 goto end of loop */
statement(4);
/* if break 3 goto end of loop */
}
/* end of loop */
statement(5);
As shown above, both types of loops can be expressed in terms of the other type of loop. However, there are certain situations where one might more appropriate to use than the other. For instance:
for loops are good for cases where there is a
known number of iterations, such as looping over the contents of
a list. while loops are good for cases where an action
should be performed so long as a condition applies, such as
reading a file while there is still data in it. do..while loops are best suited for cases
where an action needs to be performed until its results
cause a condition to be met, such as performing a calculation
until the accuracy achieved is good enough. || operator Here is my short but complete program demonstrating the
|| (logical or) operator:
main() {
int a, b;
printf("Enter two space separated numbers: ");
scanf("%d %d", &a, &b);
if (a || b) {
printf("Either %d or %d is non-zero (true).\n", a, b);
} else {
printf("Both %d and %d are zero (false).\n", a, b);
}
}
It takes two integers in the range (-2147483647-1) to (2147483647) and tells you if they are both zero of if either is non-zero. For example:
bash$ a.out Enter two space separated numbers: -2147483648 0 Either -2147483648 or 0 is non-zero (true). bash$ a.out Enter two space separated numbers: 0 0 Both 0 and 0 are zero (false). bash$
Pseudo-code:
0), give the answer
one (1). Otherwise, continue. 1). 2)
to n, multiplying result by
index at each step. Since factorials can be huuuge, special considerations will need to be made about the data type used.
/* Moderately Big Factorials; version 1.0 */
/* (c) Ian Hickson 2000, distributed under the GNU GPL */
main() {
/* Declare variables: */
unsigned long long int n, index;
long double result = 1.0;
/* ... result is floating point because it would overflow integer types */
/* Introduce program: */
printf("Welcome to the Moderately Big Factorials program, version 1.0\n");
printf("This program was written and is (c) copyright 2000 by Ian Hickson.\n");
printf("It is covered by the GNU General Public Licence.\n\n");
/* Get number of which to find factorial: */
printf("Please enter an integer in the range 0 to 1754: ");
scanf("%llu", &n);
/* ... "llu" means "expect unsigned long long int" */
/* Perform factorial by multiplying 'result' by itself
and every number from 2 to n sequentially: */
for (index = 1; index < n; result *= ++index);
/* Print output: */
printf("\n%llu! = %Lg\n", n, result);
/* ... "Lg" means "expect long double and print it nicely" */
}
Some sample output:
bash$ a.out Welcome to the Moderately Big Factorials program, version 1.0 This program was written and is (c) copyright 2000 by Ian Hickson. It is covered by the GNU General Public Licence. Please enter an integer in the range 0 to 1754: 0 0! = 1 bash$ a.out Welcome to the Moderately Big Factorials program, version 1.0 This program was written and is (c) copyright 2000 by Ian Hickson. It is covered by the GNU General Public Licence. Please enter an integer in the range 0 to 1754: 3 3! = 6 bash$ a.out Welcome to the Moderately Big Factorials program, version 1.0 This program was written and is (c) copyright 2000 by Ian Hickson. It is covered by the GNU General Public Licence. Please enter an integer in the range 0 to 1754: 1754 1754! = 1.97926e+4930 bash$
Real numbers in a computer are stored in a fixed number of
bits. For example, a number stored as a float is
stored in 32 bits, with the first bit representing the sign
(s), the next 8 representing the exponent
(e), and the last 23 representing the mantissa
(m). The number is then given by:
(-1)s × m × be
...where b is the base, namely 2. This implies that only a certain number (23 bits worth) of significant digits can be stored. So if the exponent is such that the least significant bit of the mantissa actually represents 2 (or more), then adding 1 to the number will not have any effect, and that is the limit on numerical accuracy.
These problems can be reduced somewhat by using rather large
datatypes such as the 128 bit long double.
Unfortunately, the math C library does not support
long double versions of its functions, so for instance
you cannot use exp() or pow() with
long doubles.
My integrand doesn't really take a simpler form as x gets very large. In the limit of x → ∞:
e-x2 = 0
...i.e., it is asymptotic about y=0.
/* Simple Integrator; version 1.0 */
/* (c) Ian Hickson 2000, distributed under the GNU GPL */
#include <stdio.h> /* for input/output */
#include <math.h> /* for maths */
#define MaxSignificantFigures 16
#define MinSignificantFigures 3
#define Function exp(-pow(x,2))
#define KnownAnswer 0.5*sqrt(M_PI)
#define LowLimit 0 /* lower limit -- set by brief */
#undef HighLimit /* upper limit -- use findHighLimit() */
/* cannot use "long double" mathematics because exp() and pow() don't
support long doubles!!! */
/* Wrapper around the function f(x) which we will integrate */
double f(double x) {
return Function;
}
/* Returns the number to use as the high limit, typically, the number
to use as infinity. This is very roughly the point at which the
value of the function to be integrated becomes negligible (f(xN)=0)
*/
double findHighLimit() {
#ifndef HighLimit /* otherwise we can simply use given value */
/* Sweeps values of xN to find at what point f(x) becomes negligible */
/* WARNING! This assumes that the function is asymptotic about y=0
and does not cross the x axis. More importantly, it assumes that
the high limit is infinity! #define HighLimit if these
assumptions do not apply. */
/* Note: Granularity is 1.0 */
double xN;
for (xN=LowLimit; f(xN); xN++);
return xN;
#else
return HighLimit;
#endif
}
/* Returns the number to use as the low limit, typically 0 */
double findLowLimit() {
#ifndef LowLimit
#error no low limit algorithm defined!
#else
return LowLimit;
#endif
}
/* equivalent of abs() but for doubles (abs() is for ints only) */
double absd(double x) {
return x < 0.0 ? -x : x;
}
/* returns the number of common signifcant figures between doubles a
and b */
int getCommonSignificantFigures(double a, double b) {
/* Method:
+ multiply or divide a and b until one is 0.1 <= x < 1
and the other is y < 1
+ multiply both by ten, comparing their integer portions until
they differ, or until the number of comparisons exceeds the
number of significant digits available in the type used (as
specified by MaxSignificantFigures).
Possible bugs XXX: This will overflow if a or b are close to the
maximum or minimum exponents already (won't occur with the
equation of assessment version 54).
*/
int MatchingSignificantFigures;
while (absd(a) < 0.1 || absd(b) < 0.1) {
a *= 10;
b *= 10;
}
while (absd(a) >= 1 || absd(b) >= 1) {
a /= 10;
b /= 10;
}
MatchingSignificantFigures = -1;
do {
MatchingSignificantFigures++;
a *= 10;
b *= 10;
} while (((long long int)a == (long long int)b) &&
(MatchingSignificantFigures < MaxSignificantFigures));
return MatchingSignificantFigures;
}
/* Integrate from x1 to xN in lastN steps */
/* Uses the Extended Trapezoidal Numerical method. */
double Integrate(double x1, double xN, int lastN) {
double I = 0.5 * (f(x1)+f(xN)); /* First and last terms */
double h = (xN - x1) / (lastN - 1); /* Find the trapezium width */
double stepDelta;
int i;
for (i=1; i<lastN; i++) {
stepDelta = f(x1 + i*h);
if (stepDelta + I == I) { break; } /* stop if addition is negligible */
I += stepDelta;
}
I *= h;
return I;
}
/* Perform the integration estimation by successively increasing the
number of strips used (N) until a sufficiently accurate result is
found. This is hugely innefficient. A better way of doing this
would be to only recalculate the required numbers on each pass.
However, for a function as simple as f(x) here, the overhead of
either storing the known values of f(x) or calculating which values
have already been taken care of would probably be more than the
cost of simply brute-forcing our way through. */
double EstimateIntegral(double x1, double xN, int SignificantFiguresRequired) {
double lastResult, thisResult;
int N = 2; /* Suitable initial value of N */
thisResult = Integrate(x1, xN, N); /* initialize the thisResult variable */
do {
N += 1; /* sweep through values of N linearly */
lastResult = thisResult;
thisResult = Integrate(x1, xN, N);
} while (getCommonSignificantFigures(lastResult, thisResult)
<= SignificantFiguresRequired);
/* ...note that we always require the CommonSignificantFigures to be
at least one more than the SignificantFiguresRequired. */
return thisResult;
}
/* Main program code */
void main() {
double x1, xN, Result;
int SignificantFiguresRequired;
/* Introduce program: */
printf("Welcome to the Simple Integrator program, version 1.0\n");
printf("This program was written and is (c) copyright 2000 by Ian Hickson.\n");
printf("It is covered by the GNU General Public Licence.\n\n");
/* Setup variables: */
x1 = findLowLimit();
xN = findHighLimit();
printf("Integrand: exp(-x^2)\nLimits: %g..Infinity (Infinity=%g)\n", x1, xN);
/* Find out what accuracy we need: */
printf("\nEnter the desired accuracy in terms of the number of base 10 significant\n\
figures to be confident about: ");
scanf("%d", &SignificantFiguresRequired);
/* Do it, if the range of accuracy is sensible */
if (SignificantFiguresRequired < MaxSignificantFigures) {
if (SignificantFiguresRequired > MinSignificantFigures) {
Result = EstimateIntegral(x1, xN, SignificantFiguresRequired);
printf(" It would seem the integral is: %.*g to %d significant figures\n",
SignificantFiguresRequired, Result, SignificantFiguresRequired);
printf(" My full estimate is..........: %.*g\n", MaxSignificantFigures, Result);
printf(" The actual answer is.........: %.*g to %d significant figures [KNOWN]\n",
SignificantFiguresRequired, KnownAnswer, SignificantFiguresRequired);
printf(" The full actual answer is....: %.*g [KNOWN]\n",
MaxSignificantFigures, KnownAnswer);
printf(" The first error is in significant digit number %d.\n",
getCommonSignificantFigures(KnownAnswer, Result)+1);
printf(" The error is +/-%g or %g%%.\n\n", absd(Result-KnownAnswer),
absd(100.0*((Result-KnownAnswer)/KnownAnswer)));
} else {
printf("Sorry, I cannot be confident about numbers with such a low\
accuracy (less than %d).\n", MinSignificantFigures+1);
}
} else {
printf("Sorry, I cannot be confident about that kind of high accuracy\
(more than %d).\n", MaxSignificantFigures-1);
}
}
Sample output:
bash$ a.out Welcome to the Simple Integrator program, version 1.0 This program was written and is (c) copyright 2000 by Ian Hickson. It is covered by the GNU General Public Licence. Integrand: exp(-x^2) Limits: 0..Infinity (Infinity=28) Enter the desired accuracy in terms of the number of base 10 significant figures to be confident about: 4 It would seem the integral is: 0.8862 to 4 significant figures My full estimate is..........: 0.8862313935544009 The actual answer is.........: 0.8862 to 4 significant figures [KNOWN] The full actual answer is....: 0.8862269254527579 [KNOWN] The first error is in significant digit number 5. The error is +/-4.4681e-06 or 0.0015839%. bash$ a.out Welcome to the Simple Integrator program, version 1.0 This program was written and is (c) copyright 2000 by Ian Hickson. It is covered by the GNU General Public Licence. Integrand: exp(-x^2) Limits: 0..Infinity (Infinity=28) Enter the desired accuracy in terms of the number of base 10 significant figures to be confident about: 15 It would seem the integral is: 0.886226925452758 to 15 significant figures My full estimate is..........: 0.8862269254527583 The actual answer is.........: 0.886226925452758 to 15 significant figures [KNOWN] The full actual answer is....: 0.8862269254527579 [KNOWN] The first error is in significant digit number 15. The error is +/-3.33067e-16 or 1.18069e-13%. bash$
Since I needed to do this to check the validity of the
program, the code given for questions 2c and 3 already perform the
necessary steps for finding if the program produces the estimate
to the correct level of precision. The code in question is the
getCommonSignificantFigures() function and in
particular the calling of this function in the output section of
the main() block:
printf(" The first error is in significant digit number %d.\n",
getCommonSignificantFigures(KnownAnswer, Result)+1);
See the previous answer for sample output of this code.