1 ! Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at 2 ! the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights 3 ! reserved. See files LICENSE and NOTICE for details. 4 ! 5 ! This file is part of CEED, a collection of benchmarks, miniapps, software 6 ! libraries and APIs for efficient high-order finite element and spectral 7 ! element discretizations for exascale applications. For more information and 8 ! source code availability see http://github.com/ceed. 9 ! 10 ! The CEED research is supported by the Exascale Computing Project 17-SC-20-SC, 11 ! a collaborative effort of two U.S. Department of Energy organizations (Office 12 ! of Science and the National Nuclear Security Administration) responsible for 13 ! the planning and preparation of a capable exascale ecosystem, including 14 ! software, applications, hardware, advanced system engineering and early 15 ! testbed platforms, in support of the nation's exascale computing imperative. 16 ! 17 !----------------------------------------------------------------------- 18 subroutine buildmats(qref,qweight,interp,grad) 19 integer p,q,d 20 parameter(p=6) 21 parameter(q=4) 22 parameter(d=2) 23 24 real*8 qref(d*q) 25 real*8 qweight(q) 26 real*8 interp(p*q) 27 real*8 grad(d*p*q) 28 real*8 x1,x2 29 30 qref=(/2.d-1,6.d-1,1.d0/3.d0,2.d-1,2.d-1,2.d-1,1.d0/3.d0,6.d-1/) 31 qweight=(/25.d0/96.d0,25.d0/96.d0,-27.d0/96.d0,25.d0/96.d0/) 32 33 do i=0,q-1 34 x1 = qref(0*q+i+1) 35 x2 = qref(1*q+i+1); 36 ! Interp 37 interp(i*P+1)=2.*(x1+x2-1.)*(x1+x2-1./2.); 38 interp(i*P+2)=-4.*x1*(x1+x2-1.); 39 interp(i*P+3)=2.*x1*(x1-1./2.); 40 interp(i*P+4)=-4.*x2*(x1+x2-1.); 41 interp(i*P+5)=4.*x1*x2; 42 interp(i*P+6)=2.*x2*(x2-1./2.); 43 ! Grad 44 grad((i+0)*P+1)=2.*(1.*(x1+x2-1./2.)+(x1+x2-1.)*1.); 45 grad((i+Q)*P+1)=2.*(1.*(x1+x2-1./2.)+(x1+x2-1.)*1.); 46 grad((i+0)*P+2)=-4.*(1.*(x1+x2-1.)+x1*1.); 47 grad((i+Q)*P+2)=-4.*(x1*1.); 48 grad((i+0)*P+3)=2.*(1.*(x1-1./2.)+x1*1.); 49 grad((i+Q)*P+3)=2.*0.; 50 grad((i+0)*P+4)=-4.*(x2*1.); 51 grad((i+Q)*P+4)=-4.*(1.*(x1+x2-1.)+x2*1.); 52 grad((i+0)*P+5)=4.*(1.*x2); 53 grad((i+Q)*P+5)=4.*(x1*1.); 54 grad((i+0)*P+6)=2.*0.; 55 grad((i+Q)*P+6)=2.*(1.*(x2-1./2.)+x2*1.); 56 enddo 57 58 end 59 !----------------------------------------------------------------------- 60