| 1 | #include <stdio.h>
|
|---|
| 2 | #ifdef _CIVL
|
|---|
| 3 | #include <civlc.cvh>
|
|---|
| 4 | #endif
|
|---|
| 5 |
|
|---|
| 6 | #define UPPER(up, step, offset) \
|
|---|
| 7 | up % step == offset ? up : \
|
|---|
| 8 | up % step < offset ? up + (offset - up % step) : up - (up % step - offset) + step
|
|---|
| 9 |
|
|---|
| 10 | #ifdef _CIVL
|
|---|
| 11 | $input int ni,nj;
|
|---|
| 12 | $input int bi,bj; //TODO: bounds of block sizes ?
|
|---|
| 13 | $assume(ni > 1);
|
|---|
| 14 | $assume(nj > 1);
|
|---|
| 15 | $input double uin[ni][nj];
|
|---|
| 16 | #endif
|
|---|
| 17 | double u1[ni][nj];
|
|---|
| 18 | double u2[ni][nj];
|
|---|
| 19 |
|
|---|
| 20 | int main(int argc, char** argv) {
|
|---|
| 21 | int i,j,ir,jr,ib,jb;
|
|---|
| 22 |
|
|---|
| 23 | /*
|
|---|
| 24 | for(i=0;i<ni;i++) {
|
|---|
| 25 | for(j=0;j<nj;j++) {
|
|---|
| 26 | u1[i][j] = 0;
|
|---|
| 27 | u2[i][j] = 0;
|
|---|
| 28 | }
|
|---|
| 29 | }
|
|---|
| 30 |
|
|---|
| 31 | // straightforward loop nest
|
|---|
| 32 | for(i=1;i<ni-1;i++) {
|
|---|
| 33 | for(j=1;j<nj-1;j++) {
|
|---|
| 34 | u1[i][j] = uin[i][j]*uin[i][j];
|
|---|
| 35 | }
|
|---|
| 36 | }
|
|---|
| 37 | */
|
|---|
| 38 | $havoc(&u1);
|
|---|
| 39 | $havoc(&u2);
|
|---|
| 40 | $assume($forall (int i : 0 .. ni-1) u1[i][0] == 0 && u1[i][nj-1] == 0 && u2[i][0] == 0 && u2[i][nj-1] == 0);
|
|---|
| 41 | $assume($forall (int j : 0 .. nj-1) u1[0][j] == 0 && u1[ni-1][j] == 0 && u2[0][j] == 0 && u2[ni-1][j] == 0);
|
|---|
| 42 | $assume($forall (int i : 1 .. ni-2) ($forall (int j : 1 .. nj-2) u1[i][j] == uin[i][j] * uin[i][j]));
|
|---|
| 43 |
|
|---|
| 44 | if (ni > 2) {
|
|---|
| 45 | i = ni - 1;
|
|---|
| 46 | if (nj > 2)
|
|---|
| 47 | j = nj - 1;
|
|---|
| 48 | else
|
|---|
| 49 | j = 1;
|
|---|
| 50 | } else {
|
|---|
| 51 | i = 1;
|
|---|
| 52 | j = nj;
|
|---|
| 53 | }
|
|---|
| 54 |
|
|---|
| 55 | // blocked loop nest
|
|---|
| 56 |
|
|---|
| 57 | /*@ loop invariant 1 <= i <= UPPER(ni-bi, bi, 1);
|
|---|
| 58 | @ loop invariant \forall int t; 1 <= t < i ==>
|
|---|
| 59 | @ (\forall int k; 1 <= k < UPPER(nj-bj, bj, 1) ==>
|
|---|
| 60 | @ u2[t][k] == uin[t][k]*uin[t][k]);
|
|---|
| 61 | @ loop assigns i, j, ib, jb, u2[1 .. (UPPER(ni-bi, bi, 1)-1)][1 .. (UPPER(nj-bj, bj, 1)-1)];
|
|---|
| 62 | @*/
|
|---|
| 63 | for(i=1;i<ni-bi;i=i+bi) {
|
|---|
| 64 | /*@ loop invariant 1 <= j <= UPPER(nj-bj, bj, 1);
|
|---|
| 65 | @ loop invariant \forall int t; 0 <= t < ib ==>
|
|---|
| 66 | @ (\forall int k; 1 <= k < UPPER(nj-bj, bj, 1) ==>
|
|---|
| 67 | @ u2[i+t][k] == uin[i+t][k]*uin[i+t][k]);
|
|---|
| 68 | @ loop assigns j, ib, jb, u2[i .. i + bi - 1][1 .. (UPPER(nj-bj, bj, 1)-1)];
|
|---|
| 69 | @*/
|
|---|
| 70 | for(j=1;j<nj-bj;j=j+bj) {
|
|---|
| 71 | /*@ loop invariant 0 <= ib <= bi;
|
|---|
| 72 | @ loop invariant \forall int t; 0 <= t < ib ==>
|
|---|
| 73 | @ (\forall int k; 0 <= k < bj ==> u2[i+t][j+k] == uin[i+t][j+k]*uin[i+t][j+k]);
|
|---|
| 74 | @ loop assigns ib, jb, u2[i .. i + bi - 1][j .. j + bj - 1];
|
|---|
| 75 | @*/
|
|---|
| 76 | for(ib=0;ib<bi;ib++) {
|
|---|
| 77 | /*@ loop invariant 0 <= jb <= bj;
|
|---|
| 78 | @ loop invariant \forall int k; 0 <= k < jb ==> u2[i+ib][j+k] == uin[i+ib][j+k]*uin[i+ib][j+k];
|
|---|
| 79 | @ loop assigns jb, u2[i + ib][j .. j + bj - 1];
|
|---|
| 80 | @*/
|
|---|
| 81 | for(jb=0;jb<bj;jb++) {
|
|---|
| 82 | u2[i+ib][j+jb] = uin[i+ib][j+jb]*uin[i+ib][j+jb];
|
|---|
| 83 | }
|
|---|
| 84 | }
|
|---|
| 85 | }
|
|---|
| 86 | }
|
|---|
| 87 | // remainder i
|
|---|
| 88 | for(ir=i;ir<ni-1;ir=ir++) {
|
|---|
| 89 | for(jr=1;jr<nj-bj;jr=jr++) {
|
|---|
| 90 | for(jb=0;jb<bj;jb++) {
|
|---|
| 91 | u2[ir][jr+jb] = uin[ir][jr+jb]*uin[ir][jr+jb];
|
|---|
| 92 | }
|
|---|
| 93 | }
|
|---|
| 94 | }
|
|---|
| 95 | // remainder j
|
|---|
| 96 | for(jr=j;jr<nj-1;jr=jr++) {
|
|---|
| 97 | for(ir=1;ir<ni-bi;ir=ir++) {
|
|---|
| 98 | for(ib=0;ib<bi;ib++) {
|
|---|
| 99 | u2[ir+ib][jr] = uin[ir+ib][jr]*uin[ir+ib][jr];
|
|---|
| 100 | }
|
|---|
| 101 | }
|
|---|
| 102 | }
|
|---|
| 103 | // remainder ij
|
|---|
| 104 | for(ir=i;ir<ni-1;ir=ir++) {
|
|---|
| 105 | for(jr=j;jr<nj-1;jr=jr++) {
|
|---|
| 106 | u2[ir][jr] = uin[ir][jr]*uin[ir][jr];
|
|---|
| 107 | }
|
|---|
| 108 | }
|
|---|
| 109 |
|
|---|
| 110 |
|
|---|
| 111 | // assert equality of u1 and u2
|
|---|
| 112 | for(i=0;i<ni;i++) {
|
|---|
| 113 | for(j=0;j<nj;j++) {
|
|---|
| 114 | $assert(u1[i][j] == u2[i][j]);
|
|---|
| 115 | }
|
|---|
| 116 | }
|
|---|
| 117 | }
|
|---|