/*---------------------------------------------------------------------- File : mps.c Contents: compute maximal partial sum of a sequence of numbers Author : Christian Borgelt History : 20.11.1997 file created 05.12.1997 addition counter and some comments added ----------------------------------------------------------------------*/ #include #include /*---------------------------------------------------------------------- Preprocessor Definitions ----------------------------------------------------------------------*/ #define MAXSEQLEN 10000 /* maximal length of sequence */ /*---------------------------------------------------------------------- Global Variables ----------------------------------------------------------------------*/ float seq[MAXSEQLEN]; /* buffer for number sequence */ int cnt = 0; /* counter for number of additions */ /*---------------------------------------------------------------------- Functions ----------------------------------------------------------------------*/ float mps1 (float seq[], int n) { /* --- time complexity O(n^3) */ int i, j, k; /* loop variables */ float sum; /* partial sum of sequence */ float max = 0; /* maximum of partial sums */ for (i = 0; i < n; i++) { /* traverse starting points */ for (j = i; j < n; j++) { /* traverse end points */ for (sum = 0, k = i; k <= j; k++) { sum += seq[k]; /* sum elements in section [i, j] */ cnt++; /* count the addition carried out */ } if (sum > max) max = sum; /* if the computed partial sum is */ } /* larger than the current maximum, */ } /* replace the current maximum */ return max; /* return maximal partial sum */ } /* mps1() */ /*---------------------------------------------------------------------- It is obvious that the above algorithms solves the maximal partial sum problem, since it computes all possible partial sums and determines their maximum. This is achieved by traversing all possible starting points for a sum. For each starting point all possible end points are tried. The innermost sum simply sums all the numbers between the current starting point and the current end point (including the numbers at the starting point and at the end point of the partial sequence). ----------------------------------------------------------------------*/ float mps2 (float seq[], int n) { /* --- time complexity O(n^2) */ int i, j; /* loop variables */ float sum; /* partial sum of sequence */ float max = 0; /* maximum of partial sums */ for (i = 0; i < n; i++) { /* traverse starting points */ for (sum = 0, j = i; j < n; j++) { /* traverse section */ sum += seq[j]; /* compute next partial sum */ cnt++; /* count the addition carried out */ if (sum > max) max = sum; /* if the computed partial sum is */ } /* larger than the current maximum, */ } /* replace the current maximum */ return max; /* return maximal partial sum */ } /* mps2() */ /*---------------------------------------------------------------------- The above function can be derived directly from the first version of the function `mps'. You only have to note that the innermost loop of the first version (in which k is the loop variable) computes the sums x_i, x_i +x_{i+1}, x_i +x_{i+1} +x_{i+2}, ... over and over again. The above function (second version) reuses the results of previous computations. Since the sums \sum_{k=i}^{j-1} x_k and \sum_{k=i}^j x_k differ only in the last term x_j, we need not recalculate the sum \sum_{k=i}^{j-1} x_k. We only have to add x_j to the result of the previous loop. In this way the whole loop vanishes and thus the time complexity reduces to O(n^2). ----------------------------------------------------------------------*/ float mlft (float seq[], int n) { /* --- compute maximal left sum */ int i; /* loop variable */ float sum = 0; /* current left sum */ float max = 0; /* maximal left sum */ for (i = 0; i < n; i++) { /* traverse elements from left */ sum += seq[i]; /* compute next left sum */ cnt++; /* count the addition carried out */ if (sum > max) max = sum; /* if it is greater than the current */ } /* maximum, replace the maximum */ return max; /* return maximal left sum */ } /* mlft() */ /*---------------------------------------------------------------------- `mlft' is an auxiliary function which is called from the third version of the function `mps'. It determines the maximal partial sum that starts with the first (leftmost) element of the sequence `seq'. ----------------------------------------------------------------------*/ float mrgt (float seq[], int n) { /* --- compute maximal right sum */ float sum = 0; /* current right sum */ float max = 0; /* maximal right sum */ while (--n >= 0) { /* traverse elements from right */ sum += seq[n]; /* compute next right sum */ cnt++; /* count the addition carried out */ if (sum > max) max = sum; /* if it is greater than the current */ } /* maximum, replace the maximum */ return max; /* return maximal right sum */ } /* mrgt() */ /*---------------------------------------------------------------------- `mrgt' is an auxiliary function which is called from the third version of the function `mps'. It determines the maximal partial sum that ends with the last (rightmost) element of the sequence `seq'. ----------------------------------------------------------------------*/ float mps3 (float seq[], int n) { /* --- time complexity O(n log(n)) */ int l; /* number of elements in left */ int r; /* and right part of the sequence */ float max; /* maximal partial sum */ float tmp; /* temporary buffer */ if (n <= 1) return seq[0]; /* if only one element, return it */ l = n /2; /* compute number of elements in the */ r = n -l; /* left and right part of the seq. */ max = mrgt(seq, l) /* init. max as the sum of the max. */ + mlft(seq +l, r); /* right and left sums of the parts */ cnt++; /* count the addition carried out */ tmp = mps3(seq, l); /* compute mps of left part */ if (tmp > max) max = tmp; /* and adapt the maximum */ tmp = mps3(seq +l, r); /* compute mps of right part */ if (tmp > max) max = tmp; /* and adapt the maximum */ return max; /* return maximal partial sum */ } /* mps3() */ /*---------------------------------------------------------------------- The above algorithm for the maximal partial sum problem was found by M. Shamos in 1977. It has a time complexity of O(n log(n)). See the script by Prof. Dassow for an explanation and the derivation of the time complexity. ----------------------------------------------------------------------*/ float mps4 (float seq[], int n) { /* --- time complexity O(n) */ int i; /* loop variable */ float sum = 0; /* partial sum of sequence */ float max = 0; /* maximum of partial sums */ for (i = 0; i < n; i++) { /* traverse number sequence */ sum += seq[i]; /* compute next partial sum */ cnt++; /* count the addition carried out */ if (sum < 0) /* if the preceding part of the seq. */ sum = 0; /* does not add anything, discard it */ else if (sum > max) /* if the computed partial sum is */ max = sum; /* larger than the current maximum, */ } /* replace the current maximum */ return max; /* return maximal partial sum */ } /* mps4() */ /*---------------------------------------------------------------------- The above solution for the maximal partial sum problem was found by J. Kadane in 1977. It is easy to see that it has a time complexity of O(n), but it is fairly difficult to understand why it is correct. See either the script by Prof. Dassow or the script by Prof. Ehrich for an explanation. ----------------------------------------------------------------------*/ int main (int argc, char *argv[]) { /* --- generate random numbers */ int ver; /* version of mps to run */ FILE *file = stdin; /* input file */ int n = 0; /* number of numbers read */ float res = 0; /* result of function mpsx */ if (argc != 3) { /* if wrong number of args. given */ printf("usage: %s version filename\n", argv[0]); printf("compute maximal partial sum of a sequence of numbers\n"); printf("if the file name starts with a '-', " "input is read from stdin\n"); return 0; /* print a usage message */ } /* and abort the program */ ver = atoi(argv[1]); /* get number of the version to run */ if ((ver < 1) || (ver > 4)) { /* and check this number */ printf("%s: illegal version number %d\n", argv[0], ver); return 1; } if (argv[2][0] != '-') { /* if a normal file name is given */ file = fopen(argv[2], "r"); /* try to open the input file */ if (!file) { /* if file open failed */ printf("%s: cannot open file %s\n", argv[0], argv[2]); return 1; /* print an error message */ } /* and abort the function */ } while ((n < MAXSEQLEN) /* while sequence buffer is not full */ && (fscanf(file, "%f", seq +n) == 1)) n++; /* read numbers into the vector `seq' */ if (ferror(file)) /* if a read error occurred */ printf("%s: read error on file %s\n", argv[0], argv[2]); else { /* print an error message, otherwise */ switch (ver) { /* compute maximal partial sum */ case 1: res = mps1(seq, n); break; case 2: res = mps2(seq, n); break; case 3: res = mps3(seq, n); break; case 4: res = (n <= 0) ? 0 : mps4(seq, n); break; } /* call mpsx function */ printf("%g (%d additions)\n", res, cnt); } /* print result and num. of additions */ if (file != stdin) /* if input was read from a file, */ fclose(file); /* close the file */ return 0; /* return 'ok' */ } /* main() */ /*---------------------------------------------------------------------- Note that the conditional evaluation for version 4 is necessary, since the function mps4 does not deal correctly with empty sequences (n = 0). In this program it is demonstrated how to read input from a file. Basically this is just like reading input from a terminal. Instead of `scanf' you simply use `fscanf' to read a number (the `f' indicates that this function refers to a file). `fscanf' works just like `scanf'. A format string and the address(es) of the variable(s) to store the read number(s) in have to be given. The only difference is that fscanf needs to be given a pointer to a FILE data structure as the first argument. Such a pointer can be requested by calling the `fopen' function. `fopen' simply opens a file and lets you access the file contents. It receives two arguments: the name of the file to open and a mode string. The characters in the mode string determine what operations are allowed. An `r' indicates that the file can be read, but not written, a `w' that it can written. For a more detailed explanation of the mode string see the manual page on the function `fopen'. An important thing about file operations is that there are three special files that are always open. These special files are named `stdin' (standard input), `stdout' (standard output), and `stderr' (standard error). They refer to the terminal. That is, you can read input from and write output to the terminal using file functions. For example: `fscanf(stdin, "%d", &n);' is equivalent to `scanf("%d", &n);' and `fprintf(stdout, "%d", n);' is equivalent to `printf("%d", n);'. In the above program this feature is used to make the program independent of the source of the input. It can be read from a file as well as from the terminal. This is achieved by checking the first program argument. If it does not start with a '-', it is assumed that the argument is a file name. Hence the file is opened and input is read from the file. Otherwise, that is, if the first program argument starts with a `-', input is read from the terminal. Hence no file is opened, but the special file `stdin' is used instead. After all input is read, the file has to be closed. This is done by a call to the function `fclose'. But be careful not to close one of the the special files. This can lead to serious errors. Exercise (follow up to the exercise in the program `eratos.c'): Use dynamic storage allocation to allow sequences of numbers to be processed that are longer than 10000 numbers. (Hint: The standard library function `realloc' may be useful. Start with a vector of size 0 (i.e. initialize the vector variable to NULL) and each time the vector is full and another number is read, increase the size of the vector by a reasonable amount, e.g. 1000 fields.) ----------------------------------------------------------------------*/