aboutsummaryrefslogtreecommitdiffstats
path: root/bin/rt_matrix_spin.c
diff options
context:
space:
mode:
Diffstat (limited to 'bin/rt_matrix_spin.c')
-rw-r--r--bin/rt_matrix_spin.c833
1 files changed, 0 insertions, 833 deletions
diff --git a/bin/rt_matrix_spin.c b/bin/rt_matrix_spin.c
deleted file mode 100644
index 142576c..0000000
--- a/bin/rt_matrix_spin.c
+++ /dev/null
@@ -1,833 +0,0 @@
1#include <sys/time.h>
2#include <sys/mman.h>
3
4#include <stdio.h>
5#include <stdlib.h>
6#include <unistd.h>
7#include <time.h>
8#include <string.h>
9#include <assert.h>
10#include <limits.h>
11
12
13#include "litmus.h"
14#include "common.h"
15#include "DISstressmarkRNG.h"
16
17#define MIN_SEED -2147483647
18#define MAX_SEED -1
19#define MIN_DIM 1
20#define MAX_DIM 32768
21#define MAX_ITERATIONS 65536
22#define MIN_TOLERANCE 0.000007
23#define MAX_TOLERANCE 0.5
24#define MIN_NUMBER -3.4e10/dim
25#define MAX_NUMBER 3.4e10/dim
26#define EPSI 1.0e-10
27#define MIN_DIG_NUMBER 1.0e-10
28#define MAX_DIG_NUMBER 3.4e10
29
30static char* progname;
31int loops = 1;
32//struct timeval t1, t2;
33
34/*
35 * External variable, dimension
36 */
37
38static int dim;
39
40/*
41 * matrix * vector
42 */
43
44void matrixMulvector(double *value,
45 int *col_ind,
46 int *row_start,
47 double *vector,
48 double *out)
49{
50 int l, ll;
51 double sum;
52 int tmp_rs, tmp_re;
53
54 for (l=0; l<dim; l++){
55 *(out + l) = 0;
56 tmp_rs = row_start[l];
57
58 if (tmp_rs != -1){
59 tmp_re = row_start[l+1]; /*
60 *get the start and ending elements of
61 * each row
62 */
63 for (ll=tmp_rs; ll<tmp_re; ll++){
64 *(out + l) += value[ll]*vector[col_ind[ll]];
65 }
66 }
67 }
68 return;
69}
70
71
72/*
73 * vector1 - vector2
74 */
75
76void vectorSub(double *vector1, double *vector2, double *vector){
77
78 int l;
79
80 for (l=0; l<dim; l++){
81 *(vector + l) = *(vector1 + l) - *(vector2 + l);
82 }
83 return;
84}
85
86
87/*
88 * vector1 + vector2
89 */
90
91void vectorAdd(double *vector1, double *vector2, double *vector){
92
93 int l;
94
95 for (l=0; l<dim; l++){
96 *(vector + l) = *(vector1 + l) + *(vector2 + l);
97 }
98 return;
99}
100
101/*
102 * vector1 * vector2
103 */
104
105double vectorMul(double *vector1, double *vector2){
106
107 int l;
108 double product;
109
110 product = 0;
111
112 for (l=0; l<dim; l++){
113 product += (*(vector1 + l))*(*(vector2 + l));
114
115 }
116 return product;
117}
118
119/*
120 * /vector/
121 */
122
123double vectorValue(double *vector){
124
125 double value;
126 int l;
127
128 value = 0;
129
130 for (l=0; l<dim; l++){
131 value += (*(vector + l)) * (*(vector + l));
132 }
133
134 return (sqrt(value));
135}
136
137/*
138 * transpose(vector)
139 * In fact, we return the original vector here
140 */
141
142void transpose(double *vector, double *vect){
143
144 int l;
145
146 for (l=0; l<dim; l++){
147 *(vect+l) = *(vector+l);
148 }
149 return;
150}
151
152/*
153 * value * <vector>
154 */
155void valueMulvector(double value, double *vector, double *vect){
156
157 int l;
158 int lll, i;
159 double tmp;
160
161 for (l=0; l<dim; l++){
162 *(vect + l) = *(vector + l) * value;
163 }
164 return;
165}
166
167/*
168 * generate the data distributed sparsely in matrix
169 */
170
171void initMatrix(double *matrix, int dim, int numberNonzero){
172
173 int k, l, ll;
174 int i, j;
175
176 int lll;
177 double sum;
178
179 for (k=0; k< dim*dim; k++){
180 *(matrix + k) = 0;
181 }
182
183 for (l=0; l<numberNonzero/2; l++){
184
185 i = randomUInt(1, dim-1);
186 j = randomUInt(0, i-1);
187
188 while (*(matrix + i*dim + j) != 0){
189
190 i++;
191 if (i == dim){
192 j++;
193 if (j == dim-1){
194 j = 0;
195 i = 1;
196 }
197 else{
198 i = j+1;
199 }
200 }
201 }
202
203 if (*(matrix + i*dim + j) == 0){
204 *(matrix + i*dim + j) = (double )randomNonZeroFloat(MIN_NUMBER,
205 MAX_NUMBER,
206 EPSI);
207 *(matrix + j*dim + i) = *(matrix + i*dim + j);
208 }
209 }
210
211 for (ll=0; ll<dim; ll++){
212
213
214
215 *(matrix + ll*dim + ll) = (double )randomNonZeroFloat(-MAX_DIG_NUMBER,
216 MAX_DIG_NUMBER,
217 MIN_DIG_NUMBER);
218
219 sum = 0;
220
221 for (lll=0; lll<dim; lll++){
222 if (lll != ll){
223 sum += *(matrix + lll*dim + ll);
224 }
225 }
226
227 if (*(matrix + ll*dim + ll) < sum ){
228 *(matrix + ll*dim + ll) += sum;
229 }
230 }
231
232 return;
233}
234
235/*
236 * generate the data value in the vectors
237 */
238
239void initVector(double *vector, int dim){
240
241 int l;
242
243 for (l=0; l<dim; l++){
244 *(vector + l) = (double )randomFloat (MIN_NUMBER, MAX_NUMBER);
245 }
246
247 return;
248}
249
250/*
251 * make a vector contains value of zero
252 */
253
254void zeroVector(double *vector, int dim){
255 int l;
256
257 for (l=0; l<dim; l++){
258 *(vector + l) = 0;
259 }
260 return;
261}
262
263/*
264 * return a vector which is the copy of the vect
265 */
266
267void equalVector(double *vect, double *vect1){
268
269 int l;
270
271 for (l=0; l<dim; l++){
272 *(vect1+l) = *(vect+l);
273 }
274 return;
275}
276
277
278
279void biConjugateGradient(double *value,
280 int *col_ind,
281 int *row_start,
282 double *vectorB,
283 double *vectorX,
284 double errorTolerance,
285 int maxIterations,
286 double *actualError,
287 int *actualIteration,
288 int dim)
289 /*
290 * in the code, we use a lot of temparary vectors and variables
291 * this is just for simple and clear
292 * you can optimize these temporary variables and vectors
293 * based on your need
294 *
295 */
296{
297 double *vectorR;
298 double *vectorP, *matrixAvectorP, *nextVectorR;
299 double error;
300 int iteration;
301 double alpha, beta;
302
303 double *tmpVector1, *tmpVector2, *tmpVector3;
304 double tmpValue1, tmpValue2;
305 int i;
306 int l;
307 int ll;
308
309 alpha = 0;
310 beta = 0;
311
312 vectorP = (double *)malloc(dim*sizeof(double));
313 vectorR = (double *)malloc(dim*sizeof(double));
314 nextVectorR = (double *)malloc(dim*sizeof(double));
315 vectorX = (double *)malloc(dim*sizeof(double));
316
317 tmpVector1 = (double *)malloc(dim*sizeof(double));
318 tmpVector2 = (double *)malloc(dim*sizeof(double));
319 tmpVector3 = (double *)malloc(dim*sizeof(double));
320
321 /*
322 * vectorR = vectorB - matrixA*vectorX
323 */
324 matrixMulvector(value,col_ind, row_start, vectorX, tmpVector1);
325
326 vectorSub(vectorB, tmpVector1, vectorR);
327
328 /*
329 * vectorP = vectorR
330 */
331
332 equalVector(vectorR, vectorP);
333
334 /*
335 * error = |matrixA * vectorX - vectorB| / |vectorB|
336 */
337 vectorSub(tmpVector1, vectorB, tmpVector1);
338
339 error = vectorValue(tmpVector1)/vectorValue(vectorB);
340
341 iteration = 0;
342
343 while ((iteration < maxIterations) && (error > errorTolerance)){
344
345 /*
346 * alpha = (transpose(vectorR) * vectorR) /
347 * (transpose(vectorP) * (matrixA * vectorP)
348 */
349
350 matrixMulvector(value, col_ind, row_start, vectorP, tmpVector1);
351 transpose(vectorR, tmpVector2);
352 transpose(vectorP, tmpVector3);
353 tmpValue1 = vectorMul(tmpVector3, tmpVector1);
354 tmpValue2 = vectorMul(tmpVector2, vectorR);
355 alpha = tmpValue2/tmpValue1;
356
357 /*
358 * nextVectorR = vectorR - alpha*(matrixA * vectorP)
359 */
360
361 valueMulvector(alpha, tmpVector1, tmpVector2);
362 vectorSub(vectorR, tmpVector2, tmpVector1);
363 equalVector(tmpVector1, nextVectorR);
364
365 /*
366 * beta = (transpose(nextVectorR) * nextVectorR) /
367 * (transpose(vectorR) * vectorR)
368 */
369
370 transpose(nextVectorR, tmpVector3);
371 tmpValue1 = vectorMul(tmpVector3, nextVectorR);
372 transpose(vectorR, tmpVector2);
373 tmpValue2 = vectorMul(tmpVector2, vectorR);
374 beta = tmpValue1/tmpValue2;
375
376 /*
377 * vectorX = vectorX + alpha * vectorP
378 */
379 valueMulvector(alpha, vectorP, tmpVector1);
380 vectorAdd(vectorX,tmpVector1, vectorX);
381
382 /*
383 *vectorP = nextVectorR + beta*vectorP
384 */
385 valueMulvector(beta, vectorP, tmpVector1);
386 vectorAdd(nextVectorR, tmpVector1, tmpVector1);
387
388 for (ll=0; ll<dim; ll++){
389 *(vectorP + ll) = *(tmpVector1 + ll);
390 }
391
392 /*
393 * vectorR = nextVectorR
394 */
395
396 for (l=0; l<dim; l++){
397 *(vectorR+l) = *(nextVectorR+l);
398 }
399
400 /*
401 * error = |matrixA * vectorX - vectorB| / |vectorB|
402 */
403 matrixMulvector(value, col_ind,row_start, vectorX, tmpVector1);
404 vectorSub(tmpVector1,vectorB,tmpVector1);
405 error = vectorValue(tmpVector1)/vectorValue(vectorB);
406
407 iteration++;
408 }
409
410 *actualError = error;
411 *actualIteration = iteration;
412
413 free(tmpVector1);
414 free(tmpVector2);
415 free(tmpVector3);
416
417 free(vectorR);
418 free(vectorP);
419
420 return;
421}
422
423/*
424 * This is the function to transfer the data from the matrix of dense storage
425 * to Compact Row Storage
426 */
427void create_CRS(double *matrixA,
428 double *value,
429 int *col_ind,
430 int *row_start,
431 int dim,
432 int numberNonzero)
433{
434
435 int i, j, k;
436 int cnt;
437 double tmp;
438
439 /*
440 *initialize the row_start
441 */
442
443 for(k=0; k<dim; k++){
444 row_start[k] = -1;
445 }
446
447 /*
448 * make the end of the last row to be numberNonzero + dim.
449 */
450
451 row_start[dim] = numberNonzero+dim;
452
453 /*
454 * initialize the col_ind
455 */
456
457 for (k=0; k<numberNonzero+dim; k++){
458 col_ind[k] = -1;
459 }
460
461
462 cnt = 0;
463
464 for (i=0; (cnt<numberNonzero+dim)&&(i<dim); i++){
465 for (j=0; (cnt<numberNonzero+dim)&&(j<dim); j++){
466
467 tmp = *(matrixA + i*dim + j);
468
469 if (tmp!=0){
470
471 value[cnt] = tmp;
472 col_ind[cnt] = j;
473
474 if (row_start[i] == -1)
475 row_start[i] = cnt;
476
477 cnt += 1;
478 }
479 }
480 }
481 row_start[i] = cnt;
482
483 return;
484}
485
486 int seed;
487 int numberNonzero;
488 int maxIterations;
489 float errorTolerance;
490 double actualError;
491 int actualIteration;
492
493
494 double *matrixA;
495 double *vectorB;
496 double *vectorX;
497
498 double *value;
499 int *col_ind;
500 int *row_start;
501 int sum;
502 int k;
503
504 //fscanf(stdin, "%d %d %d %d %f",
505// &seed, &dim, &numberNonzero,&maxIterations,&errorTolerance);
506int init_job() {
507
508 seed = -2;
509 //dim = 500;
510 //numberNonzero = 20000;
511 dim = randInt(100,300);
512 numberNonzero = randInt(dim+1, dim*dim/4);
513 maxIterations = 32;
514 errorTolerance = 0.002734;
515
516 assert((seed > MIN_SEED) && (seed < MAX_SEED));
517 assert((dim > MIN_DIM) && (dim < MAX_DIM));
518 assert((numberNonzero > dim) && (numberNonzero < dim*dim));
519 assert((maxIterations > 0) && (maxIterations < MAX_ITERATIONS));
520 assert((errorTolerance > MIN_TOLERANCE) && (errorTolerance < MAX_TOLERANCE));
521
522 matrixA = (double *)malloc(dim*dim*sizeof(double ));
523 vectorB = (double *)malloc(dim*sizeof(double));
524 vectorX = (double *)malloc(dim*sizeof(double));
525
526 value = (double *)malloc((numberNonzero+dim)*sizeof(double));
527 col_ind = (int *)malloc((numberNonzero+dim)*sizeof(int));
528 row_start = (int *)malloc((dim+1)*sizeof(int));
529
530 //randInit(seed);
531
532 initMatrix(matrixA, dim, numberNonzero);
533
534 create_CRS(matrixA, value, col_ind, row_start, dim, numberNonzero);
535
536 initVector(vectorB, dim);
537 zeroVector(vectorX, dim);
538
539 return 0;
540}
541
542int main_job() {
543 initVector(vectorB, dim);
544 zeroVector(vectorX, dim);
545
546 actualError = 0;
547 actualIteration = 0;
548
549 biConjugateGradient(value, col_ind, row_start, vectorB, vectorX, errorTolerance,
550 maxIterations,
551 &actualError, &actualIteration, dim);
552
553 sum = 0;
554 for (k=1; k<dim; k++){
555 sum += sum + *(vectorX + k);
556 }
557
558 //fprintf(stdout, "sum = %d, actualError = %e, actualIteration = %d\n", sum, actualError, actualIteration);
559 //fprintf(stdout, "total time = %u sec. \n", (unsigned int)endTime);
560
561 return(0);
562}
563
564int post_job() {
565 free(matrixA);
566 free(vectorB);
567 free(vectorX);
568 free(value);
569 free(col_ind);
570 free(row_start);
571 return 0;
572}
573
574static void usage(char *error) {
575 fprintf(stderr, "Error: %s\n", error);
576 fprintf(stderr,
577 "Usage:\n"
578 " rt_spin [COMMON-OPTS] WCET PERIOD DURATION\n"
579 " rt_spin [COMMON-OPTS] -f FILE [-o COLUMN] WCET PERIOD\n"
580 " rt_spin -l\n"
581 "\n"
582 "COMMON-OPTS = [-w] [-s SCALE]\n"
583 " [-p PARTITION/CLUSTER [-z CLUSTER SIZE]] [-c CLASS] [-m CRITICALITY LEVEL]\n"
584 "\n"
585 "WCET and PERIOD are microseconds, DURATION is seconds.\n");
586 exit(EXIT_FAILURE);
587}
588
589inline unsigned long get_cyclecount (void)
590{
591 unsigned long value;
592 // Read CCNT Register
593 asm volatile ("MRC p15, 0, %0, c9, c13, 0\t\n": "=r"(value));
594 return value;
595}
596
597static int loop_main(double exec_time, double emergency_exit)
598{
599 double last_loop = 0, loop_start;
600 int tmp = 0;
601
602 double start = cputime();
603 double now = cputime();
604
605 while (now + last_loop < start + exec_time) {
606 loop_start = now;
607 init_job();
608 tmp += main_job();
609 post_job();
610 now = cputime();
611 last_loop = now - loop_start;
612 if (emergency_exit && wctime() > emergency_exit) {
613 /* Oops --- this should only be possible if the execution time tracking
614 * is broken in the LITMUS^RT kernel. */
615 fprintf(stderr, "!!! rtspin/%d emergency exit!\n", getpid());
616 fprintf(stderr, "Something is seriously wrong! Do not ignore this.\n");
617 break;
618 }
619 }
620
621 return tmp;
622}
623
624static int job(double exec_time, double program_end)
625{
626 if (wctime() > program_end)
627 return 0;
628 else {
629 register int iter = 0;
630 //register unsigned long t;
631 //t = get_cyclecount();
632 //gettimeofday(&t1, NULL);
633 //while (iter++ < loops) {
634 loop_main(exec_time, program_end + 1);
635 //}
636 //t = get_cyclecount() - t;
637 //printf("%ld cycles\n", t);
638 //gettimeofday(&t2, NULL);
639 //printf("%ld us\n", ((t2.tv_sec * 1000000 + t2.tv_usec) - (t1.tv_sec * 1000000 + t1.tv_usec)));
640 sleep_next_period();
641 return 1;
642 }
643}
644
645#define OPTSTR "p:wves:l:m:i:b:"
646int main(int argc, char** argv)
647{
648 int ret;
649 lt_t wcet;
650 lt_t period;
651 lt_t budget;
652 double wcet_ms, period_ms, budget_ms;
653 unsigned int priority = LITMUS_NO_PRIORITY;
654 int migrate = 0;
655 int cluster = 0;
656 int opt;
657 int wait = 0;
658 int want_enforcement = 0;
659 double duration = 0, start = 0;
660 double scale = 1.0;
661 task_class_t class = RT_CLASS_HARD;
662 struct rt_task param;
663 struct mc2_task mc2_param;
664 struct reservation_config config;
665 int res_type = PERIODIC_POLLING;
666
667 unsigned int job_no;
668
669
670 progname = argv[0];
671
672 /* default for reservation */
673 config.id = 0;
674 config.priority = LITMUS_NO_PRIORITY; /* use EDF by default */
675 config.cpu = -1;
676
677 mc2_param.crit = CRIT_LEVEL_C;
678
679 budget_ms = 0;
680
681 while ((opt = getopt(argc, argv, OPTSTR)) != -1) {
682 switch (opt) {
683 case 'w':
684 wait = 1;
685 break;
686 case 'p':
687 cluster = atoi(optarg);
688 migrate = 1;
689 config.cpu = cluster;
690 break;
691 case 'e':
692 want_enforcement = 1;
693 break;
694 case 's':
695 scale = atof(optarg);
696 break;
697 case 'l':
698 loops = atoi(optarg);
699 break;
700 case 'm':
701 mc2_param.crit = atoi(optarg);
702 if (mc2_param.crit < CRIT_LEVEL_A || mc2_param.crit == NUM_CRIT_LEVELS) {
703 usage("Invalid criticality level.");
704 }
705 res_type = PERIODIC_POLLING;
706 break;
707 case 'b':
708 budget_ms = atof(optarg);
709 break;
710 case 'i':
711 config.priority = atoi(optarg);
712 break;
713 case ':':
714 usage("Argument missing.");
715 break;
716 case '?':
717 default:
718 usage("Bad argument.");
719 break;
720 }
721 }
722
723 if (mc2_param.crit > CRIT_LEVEL_A && config.priority != LITMUS_NO_PRIORITY)
724 usage("Bad criticailty level or priority");
725
726 if (argc - optind < 3)
727 usage("Arguments missing.");
728
729 wcet_ms = atof(argv[optind + 0]);
730 period_ms = atof(argv[optind + 1]);
731
732 wcet = ms2ns(wcet_ms);
733 period = ms2ns(period_ms);
734 if (budget_ms == 0)
735 budget_ms = wcet_ms;
736 budget = ms2ns(budget_ms);
737
738 if (wcet <= 0)
739 usage("The worst-case execution time must be a "
740 "positive number.");
741 if (period <= 0)
742 usage("The period must be a positive number.");
743 if (wcet > period) {
744 usage("The worst-case execution time must not "
745 "exceed the period.");
746 }
747
748 duration = atof(argv[optind + 2]);
749
750 if (migrate) {
751 ret = be_migrate_to_domain(cluster);
752 if (ret < 0)
753 bail_out("could not migrate to target partition or cluster.");
754 }
755
756 /* reservation config */
757 config.id = gettid();
758
759 config.polling_params.budget = budget;
760 config.polling_params.period = period;
761 config.polling_params.offset = 0;
762 config.polling_params.relative_deadline = 0;
763 if (config.polling_params.budget > config.polling_params.period) {
764 usage("The budget must not exceed the period.");
765 }
766
767 /* create a reservation */
768 ret = reservation_create(res_type, &config);
769 if (ret < 0) {
770 bail_out("failed to create reservation.");
771 }
772 //init_job();
773 srand (time(NULL));
774 randInit(-rand()%65535);
775
776 init_rt_task_param(&param);
777 param.exec_cost = wcet;
778 param.period = period;
779 param.priority = priority;
780 param.cls = class;
781 param.release_policy = TASK_PERIODIC;
782 param.budget_policy = (want_enforcement) ?
783 PRECISE_ENFORCEMENT : NO_ENFORCEMENT;
784 if (migrate) {
785 param.cpu = gettid();
786 }
787 ret = set_rt_task_param(gettid(), &param);
788
789 if (ret < 0)
790 bail_out("could not setup rt task params");
791
792 mc2_param.res_id = gettid();
793 ret = set_mc2_task_param(gettid(), &mc2_param);
794//printf("SET_MC2_TASK\n");
795 if (ret < 0)
796 bail_out("could not setup mc2 task params");
797
798 init_litmus();
799//printf("CALL\n");
800 if (mc2_param.crit == CRIT_LEVEL_C)
801 set_page_color(-1);
802 else if (mc2_param.crit < CRIT_LEVEL_C)
803 set_page_color(config.cpu);
804//printf("CALL\n");
805
806//printf("INIT_LITMUS\n");
807 start = wctime();
808 ret = task_mode(LITMUS_RT_TASK);
809//printf("TASK_MODE\n");
810 if (ret != 0)
811 bail_out("could not become RT task");
812
813
814 if (wait) {
815//printf("BEFORE WAIT\n");
816 ret = wait_for_ts_release();
817 if (ret != 0)
818 bail_out("wait_for_ts_release()");
819 start = wctime();
820 }
821
822 while (job(wcet_ms * 0.001 * scale, start + duration)) {};
823
824 ret = task_mode(BACKGROUND_TASK);
825 if (ret != 0)
826 bail_out("could not become regular task (huh?)");
827
828 reservation_destroy(gettid(), config.cpu);
829 //post_job();
830// printf("%s/%d finished.\n",progname, gettid());
831 return 0;
832}
833