1c Helper for mpn_mul_1, mpn_addmul_1, and mpn_submul_1 for Cray PVP. 2 3c Copyright 1996, 2000 Free Software Foundation, Inc. 4 5c This file is part of the GNU MP Library. 6 7c The GNU MP Library is free software; you can redistribute it and/or 8c modify it under the terms of the GNU Lesser General Public License as 9c published by the Free Software Foundation; either version 3 of the 10c License, or (at your option) any later version. 11 12c The GNU MP Library is distributed in the hope that it will be useful, 13c but WITHOUT ANY WARRANTY; without even the implied warranty of 14c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 15c Lesser General Public License for more details. 16 17c You should have received a copy of the GNU Lesser General Public License 18c along with the GNU MP Library. If not, see http://www.gnu.org/licenses/. 19 20c p1[] = hi(a[]*s); the upper limbs of each product 21c p0[] = low(a[]*s); the corresponding lower limbs 22c n is number of limbs in the vectors 23 24 subroutine gmpn_mulww(p1,p0,a,n,s) 25 integer*8 p1(0:*),p0(0:*),a(0:*),s 26 integer n 27 28 integer*8 a0,a1,a2,s0,s1,s2,c 29 integer*8 ai,t0,t1,t2,t3,t4 30 31 s0 = shiftl(and(s,4194303),24) 32 s1 = shiftl(and(shiftr(s,22),4194303),24) 33 s2 = shiftl(and(shiftr(s,44),4194303),24) 34 35 do i = 0,n-1 36 ai = a(i) 37 a0 = shiftl(and(ai,4194303),24) 38 a1 = shiftl(and(shiftr(ai,22),4194303),24) 39 a2 = shiftl(and(shiftr(ai,44),4194303),24) 40 41 t0 = i24mult(a0,s0) 42 t1 = i24mult(a0,s1)+i24mult(a1,s0) 43 t2 = i24mult(a0,s2)+i24mult(a1,s1)+i24mult(a2,s0) 44 t3 = i24mult(a1,s2)+i24mult(a2,s1) 45 t4 = i24mult(a2,s2) 46 47 p0(i)=shiftl(t2,44)+shiftl(t1,22)+t0 48 c=shiftr(shiftr(t0,22)+and(t1,4398046511103)+ 49 $ shiftl(and(t2,1048575),22),42) 50 p1(i)=shiftl(t4,24)+shiftl(t3,2)+shiftr(t2,20)+shiftr(t1,42)+c 51 end do 52 end 53