How to use libceed - 10 common examples

def test_314(ceed_resource):
  ceed = libceed.Ceed(ceed_resource)

  for dim in range (1, 4):
    P, Q = 8, 10
    Pdim = P**dim
    Qdim = Q**dim
    Xdim = 2**dim
    sum1 = sum2 = 0
    x = np.empty(Xdim*dim, dtype="float64")
    u = np.empty(Pdim, dtype="float64")

    for d in range(dim):
      for i in range(Xdim):
        x[d*Xdim + i] = 1 if (i % (2**(dim-d))) // (2**(dim-d-1)) else -1

    X = ceed.Vector(Xdim*dim)
    X.set_array(x, cmode=libceed.USE_POINTER)

nelem = 12
  dim = 2
  p = 6
  q = 4
  nx, ny = 3, 2
  ndofs = (nx*2+1)*(ny*2+1)
  nqpts = nelem*q

  # Vectors
  x = ceed.Vector(dim*ndofs)
  x_array = np.zeros(dim*ndofs)
  for i in range(ndofs):
    x_array[i] = (1. / (nx*2)) * (i % (nx*2+1))
    x_array[i+ndofs] = (1. / (ny*2)) * (i / (nx*2+1))
  x.set_array(x_array, cmode=libceed.USE_POINTER)

  qdata = ceed.Vector(nqpts)
  u = ceed.Vector(ndofs)
  v = ceed.Vector(ndofs)

  # Restrictions
  indx = np.zeros(nelem*p, dtype="int32")
  for i in range(nelem//2):
    col = i % nx;
    row = i // nx;
    offset = col*2 + row*(nx*2+1)*2

    indx[i*2*p+ 0] =  2 + offset
    indx[i*2*p+ 1] =  9 + offset
    indx[i*2*p+ 2] = 16 + offset
    indx[i*2*p+ 3] =  1 + offset

for d in range(dim):
      for i in range(Xdim):
        x[d*Xdim + i] = 1 if (i % (2**(dim-d))) // (2**(dim-d-1)) else -1

    X = ceed.Vector(Xdim*dim)
    X.set_array(x, cmode=libceed.USE_POINTER)
    Xq = ceed.Vector(Qdim*dim)
    Xq.set_value(0)
    U = ceed.Vector(Qdim)
    U.set_value(0)
    Uq = ceed.Vector(Qdim)

    bxl = ceed.BasisTensorH1Lagrange(dim, dim, 2, Q, libceed.GAUSS_LOBATTO)
    bul = ceed.BasisTensorH1Lagrange(dim, 1, Q, Q, libceed.GAUSS_LOBATTO)

    bxl.apply(1, libceed.EVAL_INTERP, X, Xq)

    xq = Xq.get_array_read()
    for i in range(Qdim):
      xx = np.empty(dim, dtype="float64")
      for d in range(dim):
        xx[d] = xq[d*Qdim + i]
      uq[i] = eval(dim, xx)

    Xq.restore_array_read()
    Uq.set_array(uq, cmode=libceed.USE_POINTER)

    # This operation is the identity because the quadrature is collocated
    bul.T.apply(1, libceed.EVAL_INTERP, Uq, U)

    bxg = ceed.BasisTensorH1Lagrange(dim, dim, 2, Q, libceed.GAUSS)
    bug = ceed.BasisTensorH1Lagrange(dim, 1, Q, Q, libceed.GAUSS)

os.path.join(file_dir, "test-qfunctions.h:setup_mass_2d"))
  qf_setup_tet.add_input("weights", 1, libceed.EVAL_WEIGHT)
  qf_setup_tet.add_input("dx", dim*dim, libceed.EVAL_GRAD)
  qf_setup_tet.add_output("rho", 1, libceed.EVAL_NONE)

  qf_mass_tet = ceed.QFunction(1, qfs.apply_mass,
                               os.path.join(file_dir, "test-qfunctions.h:apply_mass"))
  qf_mass_tet.add_input("rho", 1, libceed.EVAL_NONE)
  qf_mass_tet.add_input("u", 1, libceed.EVAL_INTERP)
  qf_mass_tet.add_output("v", 1, libceed.EVAL_INTERP)

  # Operators
  op_setup_tet = ceed.Operator(qf_setup_tet)
  op_setup_tet.set_field("weights", rxi_tet, bx_tet, libceed.VECTOR_NONE)
  op_setup_tet.set_field("dx", rx_tet, bx_tet, libceed.VECTOR_ACTIVE)
  op_setup_tet.set_field("rho", rui_tet, libceed.BASIS_COLLOCATED,
                         qdata_tet)

  op_mass_tet = ceed.Operator(qf_mass_tet)
  op_mass_tet.set_field("rho", rui_tet, libceed.BASIS_COLLOCATED, qdata_tet)
  op_mass_tet.set_field("u", ru_tet, bu_tet, libceed.VECTOR_ACTIVE)
  op_mass_tet.set_field("v", ru_tet, bu_tet, libceed.VECTOR_ACTIVE)

  ## ------------------------- Hex Elements -------------------------

  # Restrictions
  indx_hex = np.zeros(nelem_hex*p_hex*p_hex, dtype="int32")
  for i in range(nelem_hex):
    col = i % nx_hex;
    row = i // nx_hex;
    offset = (nx_tet*2+1)*(ny_tet*2)*(1+row)+col*2

import sys
from libceed import CEED_MEM_HOST, CEED_USE_POINTER
import libceed
import numpy as np

def check_values(ceed, x, value):
  b = x.get_array_read(CEED_MEM_HOST)
  for i in range(len(b)):
    if b[i] != value:
      # LCOV_EXCL_START
      print("Error reading array b[%d] = %f"%(i,b[i]))
      # LCOV_EXCL_STOP
  x.restore_array_read()

if __name__ == "__main__":
  ceed = libceed.Ceed(sys.argv[1])
  n = 10
  x = ceed.Vector(n)
  value = 1
  a = np.arange(10, 10 + n, dtype="float64")
  x.set_array(CEED_MEM_HOST, CEED_USE_POINTER, a)

  b = x.get_array_read(CEED_MEM_HOST)
  for i in range(len(b)):
    if b[i] != 10+i:
      # LCOV_EXCL_START
      print("Error reading array b[%d] = %f"%(i,b[i]))
  # LCOV_EXCL_STOP

  x.restore_array_read()

  x.set_value(3.0)

def test_505(ceed_resource):
  ceed = libceed.Ceed(ceed_resource)

  nelem = 15
  p = 5
  q = 8
  nx = nelem + 1
  nu = nelem*(p-1) + 1

  # Vectors
  x = ceed.Vector(nx)
  x_array = np.zeros(nx)
  for i in range(nx):
    x_array[i] = i / (nx - 1.0)
  x.set_array(x_array, cmode=libceed.USE_POINTER)

  qdata = ceed.Vector(nelem*q)
  u = ceed.Vector(nu)

# @file
# Test setting one vector from array of another vector

import sys
from libceed import CEED_MEM_HOST, CEED_USE_POINTER
import libceed
import numpy as np

if __name__ == "__main__":
  ceed = libceed.Ceed(sys.argv[1])

  n = 10

  x = ceed.Vector(n)
  y = ceed.Vector(n)

  a = np.arange(10, 10 + n, dtype="float64")
  x.set_array(MEM_HOST, CEED_USE_POINTER, a)

  x_array = x.GetArray(CEED_MEM_HOST)
  y.set_array(CEED_MEM_HOST, CEED_USE_POINTER, x_array)
  x.restore_array()

  y_array = y.get_array_read(CEED_MEM_HOST)
  for i in range(n):
    if y_array[i] != 10+i:

def test_304(ceed_resource, capsys):
  ceed = libceed.Ceed(ceed_resource)

  A = np.array([0.19996678, 0.0745459, -0.07448852, 0.0332866,
                0.0745459, 1., 0.16666509, -0.07448852,
                -0.07448852, 0.16666509, 1., 0.0745459,
                0.0332866, -0.07448852, 0.0745459, 0.19996678], dtype="float64")

  lam = libceed.Basis.symmetric_schur_decomposition(ceed, A, 4)

  print("Q: ")
  for i in range(4):
    for j in range(4):
      if A[j+4*i] <= 1E-14 and A[j+4*i] >= -1E-14:
         A[j+4*i] = 0
      print("%12.8f"%A[j+4*i])

  print("lambda: ")

def test_511(ceed_resource):
  ceed = libceed.Ceed(ceed_resource)

  nelem = 12
  dim = 2
  p = 6
  q = 4
  nx, ny = 3, 2
  ndofs = (nx*2+1)*(ny*2+1)
  nqpts = nelem*q

  # Vectors
  x = ceed.Vector(dim*ndofs)
  x_array = np.zeros(dim*ndofs)
  for i in range(ndofs):
    x_array[i] = (1. / (nx*2)) * (i % (nx*2+1))
    x_array[i+ndofs] = (1. / (ny*2)) * (i / (nx*2+1))
  x.set_array(x_array, cmode=libceed.USE_POINTER)

def test_500(ceed_resource):
  ceed = libceed.Ceed(ceed_resource)

  nelem = 15
  p = 5
  q = 8
  nx = nelem + 1
  nu = nelem*(p-1) + 1

  # Vectors
  x = ceed.Vector(nx)
  x_array = np.zeros(nx)
  for i in range(nx):
    x_array[i] = i / (nx - 1.0)
  x.set_array(x_array, cmode=libceed.USE_POINTER)

  qdata = ceed.Vector(nelem*q)
  u = ceed.Vector(nu)