source: CIVL/examples/pthread/CDAC/pthread-jacobi.c@ 367a390

main test-branch
Last change on this file since 367a390 was ea777aa, checked in by Alex Wilton <awilton@…>, 3 years ago

Moved examples, include, build_default.properties, common.xml, and README out from dev.civl.com into the root of the repo.

git-svn-id: ​svn://vsl.cis.udel.edu/civl/trunk@5704 fb995dde-84ed-4084-dfe6-e5aef3e2452c
  • Property mode set to 100644
File size: 8.0 KB
RevLine 
[4a64ef2]1/**************************************************************************************
2 C-DAC Tech Workshop : hyPACK-2013
3 October 15-18,2013
4 Example : pthread-jacobi.c
5
6 Objective : Jacobi method to solve AX = b matrix system of linear equations.
7
8 Input : Class Size
9 Number of Threads
10
11 Output : The solution of Ax = b or
12 The status of convergence for given bumber of iterations
13
14 Created : MAY-2013
15
16 E-mail : hpcfte@cdac.in
17
18*************************************************************************************/
19
20#include <stdio.h>
21#include<pthread.h>
22#include<sys/time.h>
23#include<stdlib.h>
[afa100f]24#include<string.h>
[4a64ef2]25
[afa100f]26#define MAX_ITERATIONS 100
[4a64ef2]27#define MAXTHREADS 8
28
29double Distance(double *X_Old, double *X_New, int matrix_size);
30pthread_mutex_t mutex1 = PTHREAD_MUTEX_INITIALIZER;
31
32double **Matrix_A, *RHS_Vector;
33double *X_New, *X_Old, *Bloc_X, rno,sum;
34
35int Number;
36void jacobi(int);
37
[afa100f]38int main(int argc, char **argv)
[4a64ef2]39{
40
41 double diag_dominant_factor = 4.0000;
42 double tolerance = 1.0E-5;
43 /* .......Variables Initialisation ...... */
44 int matrix_size, NoofRows, NoofCols,CLASS_SIZE,THREADS;
45 int NumThreads,ithread;
46 double rowsum;
47 double sum;
48 int irow, icol, index, Iteration,iteration,ret_count;
49 double time_start, time_end,memoryused;
50 struct timeval tv;
51 struct timezone tz;
52 char * CLASS;
53 FILE *fp;
54
55 pthread_attr_t pta;
56 pthread_t *threads;
57
58 memoryused =0.0;
59
60
61 printf("\n\t\t---------------------------------------------------------------------------");
62 printf("\n\t\t Centre for Development of Advanced Computing (C-DAC)");
63 printf("\n\t\t Email : hpcfte@cdac.in");
64 printf("\n\t\t---------------------------------------------------------------------------");
65 printf("\n\t\t Objective : To Solve AX=B Linear Equation (Jacobi Method)\n ");
66 printf("\n\t\t Performance for solving AX=B Linear Equation using JACOBI METHOD");
67 printf("\n\t\t..........................................................................\n");
68
69 if( argc != 3 ){
70
71 printf("\t\t Very Few Arguments\n ");
72 printf("\t\t Syntax : exec <Class-Size (Give A/B/C)> <Threads>\n");
73 exit(-1);
74 }
75 else {
[afa100f]76 CLASS = "A";
77 THREADS = 2;
[4a64ef2]78 }
79 if (THREADS > MAXTHREADS ){
80 printf("\n Number of Threads must be less than or equal to 8. Aborting ...\n");
[afa100f]81 return 0;
[4a64ef2]82 }
83 if( strcmp(CLASS, "A" )==0){
84 CLASS_SIZE = 1024;
85 }
86 if( strcmp(CLASS, "B" )==0){
87 CLASS_SIZE = 2048;
88 }
89 if( strcmp(CLASS, "C" )==0){
90 CLASS_SIZE = 4096;
91 }
92
93
94 matrix_size = CLASS_SIZE;
95 NumThreads = THREADS;
96 printf("\n\t\t Matrix Size : %d",matrix_size);
97 printf("\n\t\t Threads : %d",NumThreads);
98
99 NoofRows = matrix_size;
100 NoofCols = matrix_size;
101
102
103 /* Allocate The Memory For Matrix_A and RHS_Vector */
104 Matrix_A = (double **) malloc(matrix_size * sizeof(double *));
105 RHS_Vector = (double *) malloc(matrix_size * sizeof(double));
106
107
108 /* Populating the Matrix_A and RHS_Vector */
109 rowsum = (double) (matrix_size *(matrix_size+1)/2.0);
110 for (irow = 0; irow < matrix_size; irow++)
111 {
112 Matrix_A[irow] = (double *) malloc(matrix_size * sizeof(double));
113 sum = 0.0;
114 for (icol = 0; icol < matrix_size; icol++)
115 {
116 Matrix_A[irow][icol]= (double) (icol+1);
117 if(irow == icol ) Matrix_A[irow][icol] = (double)(rowsum);
118 sum = sum + Matrix_A[irow][icol];
119 }
120 RHS_Vector[irow] = (double)(2*rowsum) - (double)(irow+1);
121 }
122
123 memoryused+=(NoofRows * NoofCols * sizeof(double));
124 memoryused+=(NoofRows * sizeof(double));
125
126 printf("\n");
127
128 if (NoofRows != NoofCols)
129 {
130 printf("Input Matrix Should Be Square Matrix ..... \n");
131 exit(-1);
132 }
133
134 /* Dynamic Memory Allocation */
135 X_New = (double *) malloc(matrix_size * sizeof(double));
136 memoryused+=(NoofRows * sizeof(double));
137 X_Old = (double *) malloc(matrix_size * sizeof(double));
138 memoryused+=(NoofRows * sizeof(double));
139 Bloc_X = (double *) malloc(matrix_size * sizeof(double));
140 memoryused+=(NoofRows * sizeof(double));
141
142 /* Calculating the time of Operation Start */
143 gettimeofday(&tv, &tz);
144 time_start= (double)tv.tv_sec + (double)tv.tv_usec / 1000000.0;
145
146 /* Initailize X[i] = B[i] */
147
148 for (irow = 0; irow < matrix_size; irow++)
149 {
150 Bloc_X[irow] = RHS_Vector[irow];
151 X_New[irow] = RHS_Vector[irow];
152 }
153
154 /* Allocating the memory for user specified number of threads */
155 threads = (pthread_t *) malloc(sizeof(pthread_t) * NumThreads);
156
157 /* Initializating the thread attribute */
158 pthread_attr_init(&pta);
159
160 Iteration = 0;
161 do {
162 for(index = 0; index < matrix_size; index++)
163 X_Old[index] = X_New[index];
164 for(ithread=0;ithread<NumThreads;ithread++)
165 {
166
167 /* Creating The Threads */
168 ret_count=pthread_create(&threads[ithread],&pta,(void *(*) (void *))jacobi, (void *) (matrix_size/NumThreads));
169 if(ret_count)
170 {
171 printf("\n ERROR : Return code from pthread_create() is %d ",ret_count);
172 exit(-1);
173 }
174
175 }
176
177 Iteration++;
178 for (ithread=0; ithread<NumThreads; ithread++)
179 {
180 ret_count=pthread_join(threads[ithread], NULL);
181 if(ret_count)
182 {
183 printf("\n ERROR : Return code from pthread_join() is %d ",ret_count);
184 exit(-1);
185 }
186 }
187 ret_count=pthread_attr_destroy(&pta);
188 if(ret_count)
189 {
190 printf("\n ERROR : Return code from pthread_attr_destroy() is %d ",ret_count);
191 exit(-1);
192 }
193 } while ((Iteration < MAX_ITERATIONS) && (Distance(X_Old, X_New, matrix_size) >= tolerance));
194
195 /* Calculating the time at the end of Operation */
196 gettimeofday(&tv, &tz);
197 time_end= (double)tv.tv_sec + (double)tv.tv_usec / 1000000.0;
198
199 printf("\n\t\t The Jacobi Method For AX=B .........DONE");
200 printf("\n\t\t Total Number Of Iterations : %d",Iteration);
201 printf("\n\t\t Memory Utilized : %lf MB",(memoryused/(1024*1024)));
202 printf("\n\t\t Time in Seconds (T) : %lf",(time_end - time_start));
203 printf("\n\t\t..........................................................................\n");
204
205 /* Freeing Allocated Memory */
206 free(X_New);
207 free(X_Old);
208 free(Matrix_A);
209 free(RHS_Vector);
210 free(Bloc_X);
211 free(threads);
212
213}
214
215double Distance(double *X_Old, double *X_New, int matrix_size)
216{
217 int index;
218 double Sum;
219
220 Sum = 0.0;
221 for (index = 0; index < matrix_size; index++)
222 Sum += (X_New[index] - X_Old[index]) * (X_New[index] - X_Old[index]);
223 return (Sum);
224}
225
226void jacobi(int Number)
227{
228 int i,j;
229
230 for(i = 0; i < Number; i++)
231 {
232 Bloc_X[i] = RHS_Vector[i];
233
234 for (j = 0;j<i;j++)
235 {
236 Bloc_X[i] -= X_Old[j] * Matrix_A[i][j];
237 }
238 for (j = i+1;j<Number;j++)
239 {
240 Bloc_X[i] -= X_Old[j] * Matrix_A[i][j];
241 }
242 Bloc_X[i] = Bloc_X[i] / Matrix_A[i][i];
243 }
244 for(i = 0; i < Number; i++)
245
246 {
247 X_New[i] = Bloc_X[i];
248 }
249}
Note: See TracBrowser for help on using the repository browser.