<div dir="ltr"><div dir="ltr">Optimized is a configuration flag not a versions.</div><div dir="ltr"><br></div>You need to figure out your number of non-zeros per row of you global matrix, or a bound on it, and supply that in MatMPIAIJSetPreallocation. Otherwise it has to allocate and copy memory often.<div><br></div><div>You could increase your f9 on a serial run and see what runs best and then  move to parallel with a value in f6 of about 1/2 of f9.<br><div><br></div><div><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Tue, Jan 29, 2019 at 9:13 PM Yaxiong Chen <<a href="mailto:chen2018@purdue.edu" target="_blank">chen2018@purdue.edu</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">




<div dir="ltr">
<div id="gmail-m_789429084943478815gmail-m_3717110541932564007divtagdefaultwrapper" style="font-size:12pt;color:rgb(0,0,0);font-family:Calibri,Helvetica,sans-serif" dir="ltr">
<p style="margin-top:0px;margin-bottom:0px">Thanks Mark, </p>
<p style="margin-top:0px;margin-bottom:0px">    I use PETSC 3.9.4, is this the optimized version you called?</p>
<p style="margin-top:0px;margin-bottom:0px">    Actually f9 and f6 are from the PETSC example. I don't know how to set the value correctly so I commend them. The size of my elemental matrix may vary. For 2D problem, the size of elemental matrix can be 24*24 or 32*32
 or some other sizes. And the index is not continuous. In this case, the elemental matrix may interlace with each other in the global matrix, and I may have thousands of elemental matrix to be assembled. Does the preallocating suitable for this?</p>
<p style="margin-top:0px;margin-bottom:0px"><br>
</p>
<div id="gmail-m_789429084943478815gmail-m_3717110541932564007Signature">
<div id="gmail-m_789429084943478815gmail-m_3717110541932564007divtagdefaultwrapper" dir="ltr" style="font-size:12pt;color:rgb(0,0,0);background-color:rgb(255,255,255);font-family:Calibri,Arial,Helvetica,sans-serif,EmojiFont,"Apple Color Emoji","Segoe UI Emoji",NotoColorEmoji,"Segoe UI Symbol","Android Emoji",EmojiSymbols">
<p><span style="color:rgb(33,33,33);font-size:13.3333px">Yaxiong Chen, </span><br style="color:rgb(33,33,33);font-size:13.3333px">
<span style="color:rgb(33,33,33);font-size:13.3333px">Ph.D. Student </span><br style="color:rgb(33,33,33);font-size:13.3333px">
<br style="color:rgb(33,33,33);font-size:13.3333px">
<span style="font-size:13.3333px"></span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px">School of Mechanical Engineering, 3171 </span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px"><span style="color:rgb(34,34,34);font-family:Roboto,arial,sans-serif;font-size:13px">585 Purdue Mall</span></span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px">West Lafayette, IN 47907</span></p>
<div><span style="font-size:13.3333px"><br>
</span></div>
<p></p>
<p><br style="color:rgb(33,33,33);font-size:13.3333px">
<br>
</p>
</div>
</div>
<br>
<br>
<div style="color:rgb(0,0,0)">
<hr style="display:inline-block;width:98%">
<div id="gmail-m_789429084943478815gmail-m_3717110541932564007divRplyFwdMsg" dir="ltr"><font face="Calibri, sans-serif" style="font-size:11pt" color="#000000"><b>From:</b> Mark Adams <<a href="mailto:mfadams@lbl.gov" target="_blank">mfadams@lbl.gov</a>><br>
<b>Sent:</b> Tuesday, January 29, 2019 8:25 PM<br>
<b>To:</b> Yaxiong Chen<br>
<b>Cc:</b> Song Gao; <a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a><br>
<b>Subject:</b> Re: [petsc-users] PETSC matrix assembling super slow</font>
<div> </div>
</div>
<div>
<div dir="ltr">Slow assembly is often from not preallocating correctly. I am guessing that you are using Q1 element and f9==9, and thus the preallocation should be OK if this is a scalar problem on a regular grid and f6-==6 should be OK for the off processor
 allocation, if my assumptions are correct.
<div><br>
</div>
<div>You can run with -info, which will tell you how many allocation were done in assembly. Make sure that it is small (eg, 0).</div>
<div><br>
I see you use f90 array stuff 'idx-1'. Compilers can sometimes do crazy things with seeming simple code. You could just do this manually if you can find anything else.</div>
<div><br>
</div>
<div>And I trust you are running an optimized version of PETSc.<br>
<div><br>
</div>
</div>
</div>
<br>
<div class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail_quote">
<div dir="ltr" class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail_attr">On Tue, Jan 29, 2019 at 6:22 PM Yaxiong Chen via petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" id="gmail-m_789429084943478815gmail-m_3717110541932564007LPlnk396721" class="gmail-m_789429084943478815gmail-m_3717110541932564007OWAAutoLink" target="_blank">petsc-users@mcs.anl.gov</a>> wrote:<br>
</div>
<blockquote class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div dir="ltr">
<div id="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416divtagdefaultwrapper" dir="ltr" style="font-size:12pt;color:rgb(0,0,0);font-family:Calibri,Helvetica,sans-serif">
<p style="margin-top:0px;margin-bottom:0px">Hi Song,</p>
<p style="margin-top:0px;margin-bottom:0px">I don't quite understand how I can use this command. I don't partition the gloabl matrix. If I add my elemental matrix to the global system it will be like this. And in my parallel part, I use each core to generate
 the elemental matrix in turn. In this case, I guess each core will be assigned the space for global matrix and finally be assembled.But according to the manual, it seems each core will store a part of the global matrix. I<span style="font-size:12pt">s the </span><span style="font-size:12pt;font-family:Times">local
 submatrix in the </span><a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatMPIAIJSetPreallocation.html#MatMPIAIJSetPreallocation" id="gmail-m_789429084943478815gmail-m_3717110541932564007LPlnk849688" class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416OWAAutoLink gmail-m_789429084943478815gmail-m_3717110541932564007OWAAutoLink" style="font-family:Times;font-size:12pt" target="_blank">MatMPIAIJSetPreallocation</a><span style="font-family:Times;font-size:12pt">(</span><a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/Mat.html#Mat" id="gmail-m_789429084943478815gmail-m_3717110541932564007LPlnk949405" class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416OWAAutoLink gmail-m_789429084943478815gmail-m_3717110541932564007OWAAutoLink" style="font-family:Times;font-size:12pt" target="_blank">Mat</a><span style="font-family:Times;font-size:12pt">
 B,</span><a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Sys/PetscInt.html#PetscInt" id="gmail-m_789429084943478815gmail-m_3717110541932564007LPlnk752802" class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416OWAAutoLink gmail-m_789429084943478815gmail-m_3717110541932564007OWAAutoLink" style="font-family:Times;font-size:12pt" target="_blank">PetscInt</a><span style="font-family:Times;font-size:12pt">
 d_nz,const </span><a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Sys/PetscInt.html#PetscInt" id="gmail-m_789429084943478815gmail-m_3717110541932564007LPlnk373685" class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416OWAAutoLink gmail-m_789429084943478815gmail-m_3717110541932564007OWAAutoLink" style="font-family:Times;font-size:12pt" target="_blank">PetscInt</a><span style="font-family:Times;font-size:12pt">
 d_nnz[],</span><a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Sys/PetscInt.html#PetscInt" id="gmail-m_789429084943478815gmail-m_3717110541932564007LPlnk938147" class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416OWAAutoLink gmail-m_789429084943478815gmail-m_3717110541932564007OWAAutoLink" style="font-family:Times;font-size:12pt" target="_blank">PetscInt</a><span style="font-family:Times;font-size:12pt">
 o_nz,const </span><a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Sys/PetscInt.html#PetscInt" id="gmail-m_789429084943478815gmail-m_3717110541932564007LPlnk266292" class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416OWAAutoLink gmail-m_789429084943478815gmail-m_3717110541932564007OWAAutoLink" style="font-family:Times;font-size:12pt" target="_blank">PetscInt</a><span style="font-family:Times;font-size:12pt">
 o_nnz[])the same as my elemental matrix? </span></p>
<p style="margin-top:0px;margin-bottom:0px"><span style="font-family:Times"><a name="m_789429084943478815_m_3717110541932564007_x_m_9052320459455506416_MatMPIAIJSetPreallocation" id="gmail-m_789429084943478815gmail-m_3717110541932564007LPlnk139144" class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416OWAAutoLink gmail-m_789429084943478815gmail-m_3717110541932564007OWAAutoLink" style="font-family:Times"></a></span></p>
<h1><img size="3840" id="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416img659" style="max-width: 100%;" src="cid:1689caf2cd3f456b1e51"><br>
</h1>
<div><br>
</div>
<div>Thanks</div>
<p></p>
<p style="margin-top:0px;margin-bottom:0px"><br>
</p>
<div id="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416Signature">
<div id="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416divtagdefaultwrapper" dir="ltr" style="font-size:12pt;color:rgb(0,0,0);background-color:rgb(255,255,255);font-family:Calibri,Arial,Helvetica,sans-serif,EmojiFont,"Apple Color Emoji","Segoe UI Emoji",NotoColorEmoji,"Segoe UI Symbol","Android Emoji",EmojiSymbols">
<p><span style="color:rgb(33,33,33);font-size:13.3333px">Yaxiong Chen, </span><br style="color:rgb(33,33,33);font-size:13.3333px">
<span style="color:rgb(33,33,33);font-size:13.3333px">Ph.D. Student </span><br style="color:rgb(33,33,33);font-size:13.3333px">
<br style="color:rgb(33,33,33);font-size:13.3333px">
<span style="font-size:13.3333px"></span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px">School of Mechanical Engineering, 3171 </span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px"><span style="color:rgb(34,34,34);font-family:Roboto,arial,sans-serif;font-size:13px">585 Purdue Mall</span></span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px">West Lafayette, IN 47907</span></p>
<div><span style="font-size:13.3333px"><br>
</span></div>
<p></p>
<p><br style="color:rgb(33,33,33);font-size:13.3333px">
<br>
</p>
</div>
</div>
<br>
<br>
<div style="color:rgb(0,0,0)">
<hr style="display:inline-block;width:98%">
<div id="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416divRplyFwdMsg" dir="ltr"><font face="Calibri, sans-serif" color="#000000" style="font-size:11pt"><b>From:</b> Song Gao <<a href="mailto:song.gao2@mail.mcgill.ca" id="gmail-m_789429084943478815gmail-m_3717110541932564007LPlnk653762" class="gmail-m_789429084943478815gmail-m_3717110541932564007OWAAutoLink" target="_blank">song.gao2@mail.mcgill.ca</a>><br>
<b>Sent:</b> Tuesday, January 29, 2019 1:22 PM<br>
<b>To:</b> Yaxiong Chen<br>
<b>Cc:</b> Matthew Knepley; <a href="mailto:petsc-users@mcs.anl.gov" id="gmail-m_789429084943478815gmail-m_3717110541932564007LPlnk765389" class="gmail-m_789429084943478815gmail-m_3717110541932564007OWAAutoLink" target="_blank">
petsc-users@mcs.anl.gov</a><br>
<b>Subject:</b> Re: [petsc-users] PETSC matrix assembling super slow</font>
<div> </div>
</div>
<div>
<div dir="ltr">
<div dir="ltr">I think you would prefer to preallocate the matrix<br>
<br>
</div>
<div>uncomment this line <br>
</div>
<div dir="ltr">
<div>! call MatMPIAIJSetPreallocation(Amat,f9,PETSC_NULL_INTEGER,f6,PETSC_NULL_INTEGER, ierr)<br>
<br>
<br>
</div>
</div>
</div>
<br>
<div class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416x_gmail_quote">
<div dir="ltr" class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416x_gmail-m_-5873209276696769746gmail_attr">
Le mar. 29 janv. 2019, à 12 h 40, Yaxiong Chen via petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" id="gmail-m_789429084943478815gmail-m_3717110541932564007LPlnk288971" class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416OWAAutoLink gmail-m_789429084943478815gmail-m_3717110541932564007OWAAutoLink" target="_blank">petsc-users@mcs.anl.gov</a>> a écrit :<br>
</div>
<blockquote class="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416x_gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div dir="ltr">
<div id="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416x_gmail-m_-5873209276696769746gmail-m_-9202678646466270310divtagdefaultwrapper" dir="ltr" style="font-size:12pt;color:rgb(0,0,0);font-family:Calibri,Helvetica,sans-serif">
<p style="margin-top:0px;margin-bottom:0px">Hello,</p>
<p style="margin-top:0px;margin-bottom:0px"><br>
</p>
<p style="margin-top:0px;margin-bottom:0px">I have a 2D system which is assembled by each elemental matrix.  Ae is my elemental matrix, auxRHSe(:) and RHSe(:) and corresponding right hand side, idx is the global index. My code is as follow, however ,the assembling
 rate is super slow(Marked red in the code). I am not sure whether the assembling type is right or not. Since for each element, idx are not continuous numbers. Do you have any idea what is the better way to assemble the matrix?</p>
<p style="margin-top:0px;margin-bottom:0px"><br>
</p>
<p style="margin-top:0px;margin-bottom:0px">Thanks</p>
<p style="margin-top:0px;margin-bottom:0px"></p>
<p style="margin-top:0px;margin-bottom:0px"><br>
</p>
<p style="margin-top:0px;margin-bottom:0px"></p>
<div>block</div>
<div>    PetscErrorCode ierr</div>
<div>     PetscMPIInt rank,nproc, mystart</div>
<div>     PetscInt  nelem</div>
<div>     integer,allocatable  ::info(:)</div>
<div>     real(wp), allocatable :: Ae(:,:), auxRHSe(:),RHSe(:)</div>
<div>     integer, allocatable      :: idx(:)</div>
<div>     PetscScalar, pointer :: xx_v(:)</div>
<div>  </div>
<div>     PC               prec</div>
<div>     PetscScalar ::   val</div>
<div>     Vec              xvec,bvec,uvec</div>
<div>     Mat              Amat</div>
<div>     KSP              ksp</div>
<div>     PetscViewer viewer</div>
<div>     PetscInt  geq,j,k,ne,M,Istart,Iend</div>
<div>     PetscBool  flg</div>
<div>     KSPConvergedReason reason</div>
<div>     Vec         dummyVec, dummyVecs(1)</div>
<div>     MatNullSpace nullspace</div>
<div><br>
</div>
<div>     call PetscInitialize( PETSC_NULL_CHARACTER, ierr )</div>
<div><br>
</div>
<div>     if (ierr .ne. 0) then</div>
<div>     print*,'Unable to initialize PETSc'</div>
<div>     stop</div>
<div>     endif</div>
<div>     call MPI_Comm_size(PETSC_COMM_WORLD,nproc,ierr)</div>
<div>     call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)</div>
<div>     mystart=rank+1</div>
<div><br>
</div>
<div> !   Parameter set</div>
<div>     info=ptsystem%getInitInfo()</div>
<div>     nelem = info(Info_EleMatNum)+info(Info_FixedDOFNum)+info(Info_NumConstrain) </div>
<div>     print*,'nelem',nelem</div>
<div> !-------------------------------------</div>
<div> !  Create Matrix</div>
<div>     call MatCreate(PETSC_COMM_WORLD,Amat,ierr)</div>
<div>     call MatSetSizes( Amat,PETSC_DECIDE, PETSC_DECIDE, info(1), info(1), ierr )</div>
<div>     call MatSetType( Amat, MATMPIBAIJ, ierr )</div>
<div>     ! call MatMPIAIJSetPreallocation(Amat,f9,PETSC_NULL_INTEGER,f6,PETSC_NULL_INTEGER, ierr)</div>
<div>     call MatSetFromOptions( Amat, ierr )</div>
<div>     call MatSetUp( Amat, ierr )</div>
<div>     call MatGetOwnershipRange( Amat, Istart, Iend, ierr )</div>
<div> </div>
<div>      xvec = tVec(0)</div>
<div>      call MatCreateVecs( Amat, PETSC_NULL_VEC, xvec, ierr )</div>
<div>      call VecSetFromOptions( xvec, ierr )</div>
<div>      call VecDuplicate( xvec, bvec, ierr )</div>
<div>      call VecDuplicate( xvec, uvec, ierr )</div>
<div><br>
</div>
<div>   t1 = MPI_WTime(); </div>
<div><br>
</div>
<div>      do i=mystart,nelem,nproc</div>
<div>            call ptSystem%getElementalMAT(i, Ae, auxRHSe, idx)</div>
<div>            ne=size(idx)</div>
<div>            if (allocated(auxRHSe))  call VecSetValues(bvec,ne,idx-1,auxRHSe,ADD_VALUES,ierr)</div>
<div>             call MatSetValues(Amat,ne,idx-1,ne,idx-1,Ae,ADD_VALUES,ierr)</div>
<div>      end do</div>
<div><br>
</div>
<div>      nelem   = info(Info_EleRHSNum)</div>
<div>      mystart = rank+1</div>
<div><br>
</div>
<div>      do i = mystart, nelem, nproc  </div>
<div>        call ptSystem%getElementalRHS(i, RHSe, idx)</div>
<div>        print*,'idx',idx</div>
<div>        ne=size(idx)</div>
<div>        if (allocated(RHSe))  call VecSetValues(bvec,ne,idx-1,RHSe,ADD_VALUES,ierr)</div>
<div>      end do </div>
<div><span style="color:rgb(255,0,0)">      call MatAssemblyBegin(Amat,MAT_FINAL_ASSEMBLY,ierr)</span></div>
<div><span style="color:rgb(255,0,0)">      call MatAssemblyEnd(Amat,MAT_FINAL_ASSEMBLY,ierr)                          !  this part is slow, the for loop above is done but here it may get stuck</span></div>
<div><span style="color:rgb(255,0,0)">      call VecAssemblyBegin(bvec,ierr)                                                                    ! For a 2500 DOF system, assembling only takes over 2 seconds</span></div>
<div><span style="color:rgb(255,0,0)">      call VecAssemblyEnd(bvec,ierr)                                                                      ! But for a 10000 DOF system , it gets stuck</span></div>
<div>     t2 = MPI_WTime(); </div>
<div><span style="font-family:Calibri,Helvetica,sans-serif">      </span>print*,'assembeling time',t2-t1</div>
<div>      ! Solve</div>
<div>      call KSPCreate(PETSC_COMM_WORLD,ksp,ierr)</div>
<div>!  Set operators. Here the matrix that defines the linear system</div>
<div>!  also serves as the preconditioning matrix.</div>
<div>      call KSPSetOperators(ksp,Amat,Amat,ierr)</div>
<div></div>
<div>      call KSPSetFromOptions(ksp,ierr)</div>
<div>! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -</div>
<div>!                      Solve the linear system</div>
<div>! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -</div>
<div>    call KSPSetType(ksp,KSPCG,ierr)</div>
<div>    call KSPGetPC(ksp,prec,ierr)</div>
<div>   ! call KSPSetPCSide(ksp,PC_SYMMETRIC,ierr)</div>
<div>    call PCSetType(prec,PCJACOBI,ierr)</div>
<div>     call KSPSolve(ksp,bvec,xvec,ierr)</div>
<div></div>
<div>     call PetscFinalize(ierr)</div>
<div><br>
</div>
<div>end block</div>
<br>
<p></p>
<p style="margin-top:0px;margin-bottom:0px"><br>
</p>
<div id="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416x_gmail-m_-5873209276696769746gmail-m_-9202678646466270310Signature">
<div id="gmail-m_789429084943478815gmail-m_3717110541932564007x_gmail-m_9052320459455506416x_gmail-m_-5873209276696769746gmail-m_-9202678646466270310divtagdefaultwrapper" dir="ltr" style="font-size:12pt;color:rgb(0,0,0);background-color:rgb(255,255,255);font-family:Calibri,Arial,Helvetica,sans-serif,EmojiFont,"Apple Color Emoji","Segoe UI Emoji",NotoColorEmoji,"Segoe UI Symbol","Android Emoji",EmojiSymbols">
<p><span style="color:rgb(33,33,33);font-size:13.3333px">Yaxiong Chen, </span><br style="color:rgb(33,33,33);font-size:13.3333px">
<span style="color:rgb(33,33,33);font-size:13.3333px">Ph.D. Student </span><br style="color:rgb(33,33,33);font-size:13.3333px">
<br style="color:rgb(33,33,33);font-size:13.3333px">
<span style="font-size:13.3333px"></span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px">School of Mechanical Engineering, 3171 </span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px"><span style="color:rgb(34,34,34);font-family:Roboto,arial,sans-serif;font-size:13px">585 Purdue Mall</span></span></p>
<p style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:16px"><span style="font-size:13.3333px">West Lafayette, IN 47907</span></p>
<div><span style="font-size:13.3333px"><br>
</span></div>
<p></p>
<p><br style="color:rgb(33,33,33);font-size:13.3333px">
<br>
</p>
</div>
</div>
</div>
</div>
</blockquote>
</div>
</div>
</div>
</div>
</div>
</blockquote>
</div>
</div>
</div>
</div>
</div>

</blockquote></div>
</div></div></div>