1*a1db5e05SBarry Smith #define PETSCKSP_DLL 2*a1db5e05SBarry Smith 3*a1db5e05SBarry Smith /********************************bit_mask.c************************************ 4*a1db5e05SBarry Smith 5*a1db5e05SBarry Smith Author: Henry M. Tufo III 6*a1db5e05SBarry Smith 7*a1db5e05SBarry Smith e-mail: hmt@cs.brown.edu 8*a1db5e05SBarry Smith 9*a1db5e05SBarry Smith snail-mail: 10*a1db5e05SBarry Smith Division of Applied Mathematics 11*a1db5e05SBarry Smith Brown University 12*a1db5e05SBarry Smith Providence, RI 02912 13*a1db5e05SBarry Smith 14*a1db5e05SBarry Smith Last Modification: 15*a1db5e05SBarry Smith 11.21.97 16*a1db5e05SBarry Smith *********************************bit_mask.c***********************************/ 17*a1db5e05SBarry Smith #include "../src/ksp/pc/impls/tfs/tfs.h" 18*a1db5e05SBarry Smith 19*a1db5e05SBarry Smith 20*a1db5e05SBarry Smith /*********************************bit_mask.c***********************************/ 21*a1db5e05SBarry Smith PetscErrorCode bm_to_proc( char *ptr, PetscInt p_mask, PetscInt *msg_list) 22*a1db5e05SBarry Smith { 23*a1db5e05SBarry Smith PetscInt i, tmp; 24*a1db5e05SBarry Smith 25*a1db5e05SBarry Smith PetscFunctionBegin; 26*a1db5e05SBarry Smith if (msg_list) 27*a1db5e05SBarry Smith { 28*a1db5e05SBarry Smith /* low to high */ 29*a1db5e05SBarry Smith ptr+=(p_mask-1); 30*a1db5e05SBarry Smith for (i=p_mask-1;i>=0;i--) 31*a1db5e05SBarry Smith { 32*a1db5e05SBarry Smith tmp = BYTE*(p_mask-i-1); 33*a1db5e05SBarry Smith if (*ptr&BIT_0) 34*a1db5e05SBarry Smith {*msg_list = tmp; msg_list++;} 35*a1db5e05SBarry Smith if (*ptr&BIT_1) 36*a1db5e05SBarry Smith {*msg_list = tmp+1; msg_list++;} 37*a1db5e05SBarry Smith if (*ptr&BIT_2) 38*a1db5e05SBarry Smith {*msg_list = tmp+2; msg_list++;} 39*a1db5e05SBarry Smith if (*ptr&BIT_3) 40*a1db5e05SBarry Smith {*msg_list = tmp+3; msg_list++;} 41*a1db5e05SBarry Smith if (*ptr&BIT_4) 42*a1db5e05SBarry Smith {*msg_list = tmp+4; msg_list++;} 43*a1db5e05SBarry Smith if (*ptr&BIT_5) 44*a1db5e05SBarry Smith {*msg_list = tmp+5; msg_list++;} 45*a1db5e05SBarry Smith if (*ptr&BIT_6) 46*a1db5e05SBarry Smith {*msg_list = tmp+6; msg_list++;} 47*a1db5e05SBarry Smith if (*ptr&BIT_7) 48*a1db5e05SBarry Smith {*msg_list = tmp+7; msg_list++;} 49*a1db5e05SBarry Smith ptr --; 50*a1db5e05SBarry Smith } 51*a1db5e05SBarry Smith } 52*a1db5e05SBarry Smith PetscFunctionReturn(0); 53*a1db5e05SBarry Smith } 54*a1db5e05SBarry Smith 55*a1db5e05SBarry Smith /*********************************bit_mask.c***********************************/ 56*a1db5e05SBarry Smith PetscInt ct_bits( char *ptr, PetscInt n) 57*a1db5e05SBarry Smith { 58*a1db5e05SBarry Smith PetscInt i, tmp=0; 59*a1db5e05SBarry Smith 60*a1db5e05SBarry Smith PetscFunctionBegin; 61*a1db5e05SBarry Smith for(i=0;i<n;i++) 62*a1db5e05SBarry Smith { 63*a1db5e05SBarry Smith if (*ptr&128) {tmp++;} 64*a1db5e05SBarry Smith if (*ptr&64) {tmp++;} 65*a1db5e05SBarry Smith if (*ptr&32) {tmp++;} 66*a1db5e05SBarry Smith if (*ptr&16) {tmp++;} 67*a1db5e05SBarry Smith if (*ptr&8) {tmp++;} 68*a1db5e05SBarry Smith if (*ptr&4) {tmp++;} 69*a1db5e05SBarry Smith if (*ptr&2) {tmp++;} 70*a1db5e05SBarry Smith if (*ptr&1) {tmp++;} 71*a1db5e05SBarry Smith ptr++; 72*a1db5e05SBarry Smith } 73*a1db5e05SBarry Smith 74*a1db5e05SBarry Smith return(tmp); 75*a1db5e05SBarry Smith } 76*a1db5e05SBarry Smith 77*a1db5e05SBarry Smith /*********************************bit_mask.c***********************************/ 78*a1db5e05SBarry Smith PetscInt 79*a1db5e05SBarry Smith div_ceil( PetscInt numer, PetscInt denom) 80*a1db5e05SBarry Smith { 81*a1db5e05SBarry Smith PetscInt rt_val; 82*a1db5e05SBarry Smith 83*a1db5e05SBarry Smith if ((numer<0)||(denom<=0)) 84*a1db5e05SBarry Smith {SETERRQ2(PETSC_ERR_PLIB,"div_ceil() :: numer=%D ! >=0, denom=%D ! >0",numer,denom);} 85*a1db5e05SBarry Smith 86*a1db5e05SBarry Smith /* if integer division remainder then increment */ 87*a1db5e05SBarry Smith rt_val = numer/denom; 88*a1db5e05SBarry Smith if (numer%denom) 89*a1db5e05SBarry Smith {rt_val++;} 90*a1db5e05SBarry Smith 91*a1db5e05SBarry Smith return(rt_val); 92*a1db5e05SBarry Smith } 93*a1db5e05SBarry Smith 94*a1db5e05SBarry Smith /*********************************bit_mask.c***********************************/ 95*a1db5e05SBarry Smith PetscInt 96*a1db5e05SBarry Smith len_bit_mask( PetscInt num_items) 97*a1db5e05SBarry Smith { 98*a1db5e05SBarry Smith PetscInt rt_val, tmp; 99*a1db5e05SBarry Smith 100*a1db5e05SBarry Smith if (num_items<0) 101*a1db5e05SBarry Smith {SETERRQ(PETSC_ERR_PLIB,"Value Sent To len_bit_mask() Must be >= 0!");} 102*a1db5e05SBarry Smith 103*a1db5e05SBarry Smith /* mod BYTE ceiling function */ 104*a1db5e05SBarry Smith rt_val = num_items/BYTE; 105*a1db5e05SBarry Smith if (num_items%BYTE) 106*a1db5e05SBarry Smith {rt_val++;} 107*a1db5e05SBarry Smith 108*a1db5e05SBarry Smith /* make mults of sizeof int */ 109*a1db5e05SBarry Smith if ((tmp=rt_val%sizeof(PetscInt))) 110*a1db5e05SBarry Smith {rt_val+=(sizeof(PetscInt)-tmp);} 111*a1db5e05SBarry Smith 112*a1db5e05SBarry Smith return(rt_val); 113*a1db5e05SBarry Smith } 114*a1db5e05SBarry Smith 115*a1db5e05SBarry Smith /*********************************bit_mask.c***********************************/ 116*a1db5e05SBarry Smith PetscErrorCode set_bit_mask( PetscInt *bm, PetscInt len, PetscInt val) 117*a1db5e05SBarry Smith { 118*a1db5e05SBarry Smith PetscInt i, offset; 119*a1db5e05SBarry Smith char mask = 1; 120*a1db5e05SBarry Smith char *cptr; 121*a1db5e05SBarry Smith 122*a1db5e05SBarry Smith 123*a1db5e05SBarry Smith if (len_bit_mask(val)>len) 124*a1db5e05SBarry Smith {SETERRQ(PETSC_ERR_PLIB,"The Bit Mask Isn't That Large!");} 125*a1db5e05SBarry Smith 126*a1db5e05SBarry Smith cptr = (char *) bm; 127*a1db5e05SBarry Smith 128*a1db5e05SBarry Smith offset = len/sizeof(PetscInt); 129*a1db5e05SBarry Smith for (i=0;i<offset;i++) 130*a1db5e05SBarry Smith {*bm=0; bm++;} 131*a1db5e05SBarry Smith 132*a1db5e05SBarry Smith offset = val%BYTE; 133*a1db5e05SBarry Smith for (i=0;i<offset;i++) 134*a1db5e05SBarry Smith {mask <<= 1;} 135*a1db5e05SBarry Smith 136*a1db5e05SBarry Smith offset = len - val/BYTE - 1; 137*a1db5e05SBarry Smith cptr[offset] = mask; 138*a1db5e05SBarry Smith PetscFunctionReturn(0); 139*a1db5e05SBarry Smith } 140*a1db5e05SBarry Smith 141*a1db5e05SBarry Smith /*********************************bit_mask.c***********************************/ 142*a1db5e05SBarry Smith PetscInt len_buf(PetscInt item_size, PetscInt num_items) 143*a1db5e05SBarry Smith { 144*a1db5e05SBarry Smith PetscInt rt_val, tmp; 145*a1db5e05SBarry Smith 146*a1db5e05SBarry Smith PetscFunctionBegin; 147*a1db5e05SBarry Smith rt_val = item_size * num_items; 148*a1db5e05SBarry Smith 149*a1db5e05SBarry Smith /* double precision align for now ... consider page later */ 150*a1db5e05SBarry Smith if ((tmp = (rt_val%(PetscInt)sizeof(double)))) 151*a1db5e05SBarry Smith {rt_val += (sizeof(double) - tmp);} 152*a1db5e05SBarry Smith 153*a1db5e05SBarry Smith return(rt_val); 154*a1db5e05SBarry Smith } 155*a1db5e05SBarry Smith 156*a1db5e05SBarry Smith 157*a1db5e05SBarry Smith 158