1
"""Code generation for evaluation of derivatives of finite element basis values.
2
This module generates code which is more or less a C++ representation of the code
5
__author__ = "Kristian B. Oelgaard (k.b.oelgaard@gmail.com)"
6
__date__ = "2007-04-16"
7
__copyright__ = "Copyright (C) 2007-2010 Kristian B. Oelgaard"
8
__license__ = "GNU GPL version 3 or any later version"
10
# Last changed: 2010-01-30
17
from ffc.log import error, ffc_assert
18
from ffc.evaluatebasis import _map_dof
19
from ffc.evaluatebasis import _compute_basisvalues
20
from ffc.evaluatebasis import _tabulate_coefficients
21
from ffc.cpp import tabulate_matrix, IndentControl, remove_unused, tabulate_vector
22
from ffc.quadrature.quadraturegenerator_utils import generate_loop
23
#from ffc.quadrature.symbolics import set_format
26
#from cpp import format_old as format
27
from ffc.cpp import format
29
def _evaluate_basis_derivatives_all(data_list):
30
"""Like evaluate_basis, but return the values of all basis functions (dofs)."""
32
if isinstance(data_list, str):
33
return format["exception"]("evaluate_basis_derivatives_all: %s" % data_list)
35
f_r, f_s = format["free indices"][:2]
36
f_assign = format["assign"]
39
Indent = IndentControl()
42
# FIXME: KBO: Figure out which return format to use, either:
43
# [dN0[0]/dx, dN0[0]/dy, dN0[1]/dx, dN0[1]/dy, dN1[0]/dx, dN1[0]/dy, dN1[1]/dx, dN1[1]/dy, ...]
45
# [dN0[0]/dx, dN1[0]/dx, ..., dN0[1]/dx, dN1[1]/dx, ..., dN0[0]/dy, dN1[0]/dy, ..., dN0[1]/dy, dN1[1]/dy, ...]
47
# [dN0[0]/dx, dN0[1]/dx, ..., dN1[0]/dx, dN1[1]/dx, ..., dN0[0]/dy, dN0[1]/dy, ..., dN1[0]/dy, dN1[1]/dy, ...]
48
# for vector (tensor elements), currently returning option 1.
50
# FIXME: KBO: For now, just call evaluate_basis_derivatives and map values
51
# accordingly, this will keep the amount of code at a minimum. If it turns
52
# out that speed is an issue (overhead from calling evaluate_basis), we can
53
# easily generate all the code.
55
# Get total value shape and space dimension for entire element (possibly mixed).
56
value_shape = sum(sum(data["value_shape"] or (1,)) for data in data_list)
57
space_dimension = sum(data["space_dimension"] for data in data_list)
59
# Special case where space dimension is one (constant elements)
60
if space_dimension == 1:
61
code += [format["comment"]("Element is constant, calling evaluate_basis_derivatives.")]
62
code += ["evaluate_basis_derivatives(0, n, %s, coordinates, c);" % format["argument values"]]
63
return "\n".join(code)
65
# Compute number of derivatives
66
# Get the topological dimension.
67
topological_dimension = data_list[0]["topological_dimension"]
68
code += _compute_num_derivatives(topological_dimension, Indent, format)
70
# Declare helper value to hold single dof values and reset
71
code += ["", format["comment"]("Helper variable to hold values of a single dof.")]
72
if (value_shape == 1):
73
num_vals = format["num derivatives"]
75
# FIXME: KBO: Should the str(int()) be in format?
76
num_vals = format["multiply"]([str(int(value_shape)), format["num derivatives"]])
77
code += [f_assign(Indent.indent(format["float declaration"] + format["pointer"] + "dof_values"),\
78
format["component"](format["new"] + format["float declaration"], num_vals))]
79
loop_vars = [(f_r, 0, num_vals)]
80
line = [f_assign(format["component"]("dof_values", f_r), format["floating point"](0.0))]
81
code += generate_loop(line, loop_vars, Indent, format)
83
# Create loop over dofs that calls evaluate_basis_derivatives for a single dof and
84
# inserts the values into the global array.
85
code += ["", format["comment"]("Loop dofs and call evaluate_basis_derivatives.")]
87
loop_vars_r = [(f_r, 0, space_dimension)]
89
# FIXME: KBO: Move evaluate_basis_derivatives string to cpp.py
90
lines_r += ["evaluate_basis_derivatives(%s, n, dof_values, coordinates, c);" % f_r]
92
loop_vars_s = [(f_s, 0, num_vals)]
93
index = format["add"]([format["multiply"]([f_r, num_vals]), f_s])
94
name = format["component"](format["argument values"], index)
95
value = format["component"]("dof_values", f_s)
96
lines_s = [f_assign(name, value)]
97
lines_r += generate_loop(lines_s, loop_vars_s, Indent, format)
98
code += generate_loop(lines_r, loop_vars_r, Indent, format)
100
code += ["", format["comment"]("Delete pointer.")]
101
code += [Indent.indent(format["delete pointer"]("dof_values", ""))]
103
# Generate bode (no need to remove unused)
104
return "\n".join(code)
106
def _evaluate_basis_derivatives(data_list):
107
"""Evaluate the derivatives of an element basisfunction at a point. The values are
108
computed as in FIAT as the dot product of the coefficients (computed at compile time)
109
and basisvalues which are dependent on the coordinate and thus have to be computed at
112
Currently the following elements are supported in 1D:
114
Lagrange + mixed/vector valued
115
Discontinuous Lagrange + mixed/vector valued
117
Currently the following elements are supported in 2D and 3D:
119
Lagrange + mixed/vector valued
120
Discontinuous Lagrange + mixed/vector valued
121
Crouzeix-Raviart + mixed/vector valued
122
Brezzi-Douglas-Marini + mixed/vector valued
124
Not supported in 2D or 3D:
125
Raviart-Thomas ? (not tested since it is broken in FFC, but should work)
128
if isinstance(data_list, str):
129
return format["exception"]("evaluate_basis_derivatives: %s" % data_list)
131
# Init return code and indent object
133
Indent = IndentControl()
135
# Get the element cell domain, geometric and topological dimension.
136
element_cell_domain = data_list[0]["cell_domain"]
137
geometric_dimension = data_list[0]["geometric_dimension"]
138
topological_dimension = data_list[0]["topological_dimension"]
140
# Get code snippets for Jacobian, Inverse of Jacobian and mapping of
141
# coordinates from physical element to the FIAT reference element.
142
# FIXME: KBO: Change this when supporting R^2 in R^3 elements.
143
code += [Indent.indent(format["jacobian and inverse"](geometric_dimension))]
144
code += ["", Indent.indent(format["fiat coordinate map"](element_cell_domain))]
146
# Compute number of derivatives that has to be computed, and declare an array to hold
147
# the values of the derivatives on the reference element
149
code += _compute_num_derivatives(topological_dimension, Indent, format)
151
# Generate all possible combinations of derivatives
152
code += _generate_combinations(topological_dimension, Indent, format)
154
# Generate the transformation matrix
155
code += _generate_transform(element_cell_domain, Indent, format)
158
code += _reset_values(data_list, Indent, format)
160
if len(data_list) == 1:
163
# Map degree of freedom to local degree.
164
code += [_map_dof(0, Indent, format)]
166
# Generate element code.
167
code += _generate_element_code(data, 0, Indent, format)
169
# If the element is of type MixedElement (including Vector- and TensorElement).
171
# Generate element code, for all sub-elements.
172
code += _mixed_elements(data_list, Indent, format)
173
code = remove_unused("\n".join(code))
176
def _compute_num_derivatives(topological_dimension, Indent, format):
177
"Computes the number of derivatives of order 'n' as: element.cell_shape()^n."
179
code = [format["comment"]("Compute number of derivatives.")]
180
# FIXME: KBO: Should the str(int()) be in format?
181
code.append(format["declaration"](format["uint declaration"], format["num derivatives"], str(int(1))))
183
loop_vars = [(format["free indices"][0], 0, format["argument derivative order"])]
184
# FIXME: KBO: Should the str(int()) be in format?
185
lines = [format["imul"](format["num derivatives"], str(int(topological_dimension)))]
187
code += generate_loop(lines, loop_vars, Indent, format)
191
def _generate_combinations(topological_dimension, Indent, format):
192
"Generate all possible combinations of derivatives of order 'n'"
194
# Use code from format
195
code = ["", Indent.indent(format["combinations"]\
196
% {"combinations": format["derivative combinations"],\
197
"topological_dimension-1": topological_dimension-1,\
198
"num_derivatives" : format["num derivatives"],\
199
"n": format["argument derivative order"]})]
202
def _generate_transform(element_cell_domain, Indent, format):
203
"""Generate the transformation matrix, whic is used to transform derivatives from reference
204
element back to the physical element."""
206
# Generate code to construct the inverse of the Jacobian
207
if (element_cell_domain in ["interval", "triangle", "tetrahedron"]):
208
code = ["", Indent.indent(format["transform snippet"][element_cell_domain]\
209
% {"transform": format["transform matrix"],\
210
"num_derivatives" : format["num derivatives"],\
211
"n": format["argument derivative order"],\
212
"combinations": format["derivative combinations"],\
213
"K":format["transform Jinv"]})]
215
error("Cannot generate transform for shape: %s" % element_cell_domain)
219
def _reset_values(data_list, Indent, format):
220
"Reset all components of the 'values' array as it is a pointer to an array."
221
f_assign = format["assign"]
222
code = ["", Indent.indent(format["comment"]("Reset values. Assuming that values is always an array."))]
224
# Get value shape and reset values. This should also work for TensorElement,
225
# scalar are empty tuples, therefore (1,) in which case value_shape = 1.
226
value_shape = sum(sum(data["value_shape"] or (1,)) for data in data_list)
228
# Only multiply by value shape if different from 1.
230
num_vals = format["num derivatives"]
232
# FIXME: KBO: Should the str(int()) be in format?
233
num_vals = format["multiply"]([str(int(value_shape)), format["num derivatives"]])
234
name = format["component"](format["argument values"], format["free indices"][0])
235
loop_vars = [(format["free indices"][0], 0, num_vals)]
236
lines = [f_assign(name, format["floating point"](0))]
237
code += generate_loop(lines, loop_vars, Indent, format)
241
def _generate_element_code(data, sum_value_dim, Indent, format):
242
"Generate code for each basis element"
246
# Compute basisvalues, from evaluatebasis.py
247
code += _compute_basisvalues(data, Indent, format)
249
# Tabulate coefficients
250
code += _tabulate_coefficients(data, Indent, format)
252
# Tabulate coefficients for derivatives
253
code += _tabulate_dmats(data, Indent, format)
255
# Compute the derivatives of the basisfunctions on the reference (FIAT) element,
256
# as the dot product of the new coefficients and basisvalues
257
code += _compute_reference_derivatives(data, Indent, format)
259
# Transform derivatives to physical element by multiplication with the transformation matrix
260
code += _transform_derivatives(data, sum_value_dim, Indent, format)
263
code += _delete_pointers(data, Indent, format)
267
def _mixed_elements(data_list, Indent, format):
268
"Generate code for each sub-element in the event of vector valued elements or mixed elements"
270
# Prefetch formats to speed up code generation
271
f_dof_map_if = format["dof map if"]
280
# Generate code for each element
281
for data in data_list:
283
# Get value and space dimension (should be tensor ready).
284
value_dim = sum(data["value_shape"] or (1,))
285
space_dim = data["space_dimension"]
287
# Generate map from global to local dof
288
element_code = [_map_dof(sum_space_dim, Indent, format)]
290
# Generate code for basis element
291
element_code += _generate_element_code(data, sum_value_dim, Indent, format)
293
# Increase indentation, indent code and decrease indentation.
295
if_code = remove_unused(Indent.indent("\n".join(element_code)))
296
# if_code = Indent.indent("\n".join(element_code))
299
# Create if statement and add to code.
300
code += [f_if(f_dof_map_if(sum_space_dim, sum_space_dim + space_dim -1), if_code)]
302
# Increase sum of value dimension, and space dimension
303
sum_value_dim += value_dim
304
sum_space_dim += space_dim
308
def _tabulate_dmats(data, Indent, format):
309
"Tabulate the derivatives of the polynomial base"
313
# Prefetch formats to speed up code generation
314
f_table = format["static const float declaration"]
315
f_dmats = format["dmats"]
316
f_component = format["component"]
317
f_assign = format["assign"]
319
# Get derivative matrices (coefficients) of basis functions, computed by FIAT at compile time
320
derivative_matrices = data["dmats"]
322
code += [Indent.indent(format["comment"]("Tables of derivatives of the polynomial base (transpose)."))]
324
# Generate tables for each spatial direction
325
for i, dmat in enumerate(derivative_matrices):
327
# Extract derivatives for current direction (take transpose, FIAT_NEW PolynomialSet.tabulate())
328
matrix = numpy.transpose(dmat)
330
# Get shape and check dimension (This is probably not needed)
331
shape = numpy.shape(matrix)
332
ffc_assert(shape[0] == shape[1] == data["num_expansion_members"], "Something is wrong with the shape of dmats.")
334
# Declare varable name for coefficients
335
name = f_component(f_table + f_dmats(i), [shape[0], shape[1]])
336
value = tabulate_matrix(matrix, format)
337
code += [f_assign(Indent.indent(name), Indent.indent(value)), ""]
341
def _reset_dmats(shape_dmats, indices, Indent, format):
342
f_assign = format["assign"]
343
code = [format["comment"]("Resetting dmats values to compute next derivative.")]
345
loop_vars = [(indices[0], 0, shape_dmats[0]), (indices[1], 0, shape_dmats[1])]
346
dmats_old = format["component"](format["dmats"](""), [indices[0], indices[1]])
347
lines = [f_assign(dmats_old, format["floating point"](0.0))]
348
lines += [format["if"](indices[0] + format["is equal"] + indices[1],\
349
format["assign"](Indent.indent(dmats_old), format["floating point"](1.0)))]
350
code += generate_loop(lines, loop_vars, Indent, format)
353
def _update_dmats(shape_dmats, indices, Indent, format):
354
f_assign = format["assign"]
355
code = [format["comment"]("Updating dmats_old with new values and resetting dmats.")]
356
dmats = format["component"](format["dmats"](""), [indices[0], indices[1]])
357
dmats_old = format["component"](format["dmats old"], [indices[0], indices[1]])
358
loop_vars = [(indices[0], 0, shape_dmats[0]), (indices[1], 0, shape_dmats[1])]
359
lines = [f_assign(dmats_old, dmats), f_assign(dmats, format["floating point"](0.0))]
360
code += generate_loop(lines, loop_vars, Indent, format)
363
def _compute_dmats(num_dmats, shape_dmats, available_indices, deriv_index, Indent, format):
365
s, t, u = available_indices
368
code = _reset_dmats(shape_dmats, [t, u], Indent, format)
369
code += ["", format["comment"]("Looping derivative order to generate dmats.")]
371
# Set dmats matrix equal to dmats_old
372
lines = _update_dmats(shape_dmats, [t, u], Indent, format)
373
loop_vars = [(s, 0, format["argument derivative order"])]
375
lines += ["", format["comment"]("Update dmats using an inner product.")]
376
# Create dmats matrix by multiplication
377
comb = format["component"](format["derivative combinations"], [deriv_index, s])
378
for i in range(num_dmats):
379
lines += _dmats_product(shape_dmats, comb, i, [t, u], Indent, format)
381
code += generate_loop(lines, loop_vars, Indent, format)
385
def _dmats_product(shape_dmats, index, i, indices, Indent, format):
389
loop_vars = [(t, 0, shape_dmats[0]), (u, 0, shape_dmats[1])]
390
dmats = format["component"](format["dmats"](""), [t, u])
391
dmats_old = format["component"](format["dmats old"], [tu, u])
392
value = format["multiply"]([format["component"](format["dmats"](i), [t, tu]), dmats_old])
393
name = Indent.indent(format["iadd"](dmats, value))
394
lines = generate_loop([name], [(tu, 0, shape_dmats[0])], Indent, format)
395
code = [format["if"](index + format["is equal"] + str(i),\
396
"\n".join(generate_loop(lines, loop_vars, Indent, format)))]
400
def _matrix_index(i, j, range_j):
401
"Compute the index in an array from indices in a matrix i.e., m[i][j] -> a[i*range(j)+j]"
405
access = format["add"]([range_j, j])
407
# FIXME: KBO: Should the str(int()) be in format?
408
irj = format["mul"]([str(i), range_j])
409
access = format["add"]([irj, j])
412
def _compute_reference_derivatives(data, Indent, format):
413
"""Compute derivatives on the reference element by recursively multiply coefficients with
414
the relevant derivatives of the polynomial base until the requested order of derivatives
415
has been reached. After this take the dot product with the basisvalues."""
419
# Prefetch formats to speed up code generation
420
f_comment = format["comment"]
421
f_float = format["float declaration"]
422
f_component = format["component"]
423
f_tmp = format["tmp ref value"]
424
f_dmats = format["dmats"]
425
f_assign = format["assign"]
426
f_iadd = format["iadd"]
428
f_r, f_s, f_t, f_u = format["free indices"]
430
# Get number of components, change for tensor valued elements.
431
shape = data["value_shape"]
434
elif len(shape) == 1:
435
num_components = shape[0]
437
error("Tensor valued elements are not supported yet: %s" % data["family"])
439
# Get shape of derivative matrix (they should all have the same shape).
440
shape_dmats = numpy.shape(data["dmats"][0])
442
code += [Indent.indent(f_comment("Compute reference derivatives"))]
444
# Declare pointer to array that holds derivatives on the FIAT element
445
code += [Indent.indent(f_comment("Declare pointer to array of derivatives on FIAT element"))]
446
# The size of the array of reference derivatives is equal to the number of derivatives
447
# times the number of components of the basis element
448
if (num_components == 1):
449
num_vals = format["num derivatives"]
451
# FIXME: KBO: Should the str(int()) be in format?
452
num_vals = format["multiply"]([str(int(num_components)), format["num derivatives"]])
453
code += [f_assign(Indent.indent(f_float + format["pointer"] + format["reference derivatives"]),\
454
format["component"](format["new"] + f_float, num_vals))]
455
# Reset values of reference derivatives.
456
name = format["component"](format["reference derivatives"], f_r)
457
lines = [f_assign(name, format["floating point"](0))]
458
code += generate_loop(lines, [(f_r, 0, num_vals)], Indent, format)
462
# Declare matrix of dmats (which will hold the matrix product of all combinations)
463
# and dmats_old which is needed in order to perform the matrix product.
464
code += [Indent.indent(f_comment("Declare derivative matrix (of polynomial basis)."))]
465
matrix = numpy.eye(shape_dmats[0])
466
name = f_component(f_float + f_dmats(""), [shape_dmats[0], shape_dmats[1]])
467
value = tabulate_matrix(matrix, format)
468
code += [f_assign(Indent.indent(name), Indent.indent(value)), ""]
469
code += [Indent.indent(f_comment("Declare (auxiliary) derivative matrix (of polynomial basis)."))]
470
name = f_component(f_float + format["dmats old"], [shape_dmats[0], shape_dmats[1]])
471
code += [f_assign(Indent.indent(name), Indent.indent(value)), ""]
473
# Loop all derivatives and compute value of the derivative as:
474
# deriv_on_ref[r] = coeff[dof][s]*dmat[s][t]*basis[t]
475
code += [Indent.indent(f_comment("Loop possible derivatives."))]
476
loop_vars = [(f_r, 0, format["num derivatives"])]
478
# Compute dmats as a recursive matrix product
479
lines = _compute_dmats(len(data["dmats"]), shape_dmats, [f_s, f_t, f_u], f_r, Indent, format)
481
# Compute derivatives for all components
482
loop_vars_c = [(f_s, 0, shape_dmats[0]),(f_t, 0, shape_dmats[1])]
484
for i in range(num_components):
485
comp_access = _matrix_index(i, f_r, format["num derivatives"])
486
name = format["component"](format["reference derivatives"], comp_access)
487
coeffs = format["component"](format["coefficients"](i), [format["local dof"], f_s])
488
dmats = format["component"](format["dmats"](""), [f_s, f_t])
489
basis = format["component"](format["basisvalues"], f_t)
490
value = format["multiply"]([coeffs, dmats, basis])
491
lines_c.append(f_iadd(name, value))
492
lines += generate_loop(lines_c, loop_vars_c, Indent, format)
494
# Apply transformation if applicable.
495
mapping = data["mapping"]
496
if mapping == "affine":
498
elif mapping == "contravariant piola":
499
lines += ["", Indent.indent(format["comment"]\
500
("Using contravariant Piola transform to map values back to the physical element"))]
501
# Get temporary values before mapping.
502
lines += [format["const float declaration"](Indent.indent(f_tmp(i)),\
503
f_component(format["reference derivatives"], _matrix_index(i, f_r, format["num derivatives"]))) for i in range(num_components)]
505
# Create names for inner product.
506
topological_dimension = data["topological_dimension"]
507
basis_col = [f_tmp(j) for j in range(topological_dimension)]
508
for i in range(num_components):
510
jacobian_row = [format["transform"]("J", i, j, None) for j in range(topological_dimension)]
512
# Create inner product and multiply by inverse of Jacobian.
513
inner = [format["multiply"]([jacobian_row[j], basis_col[j]]) for j in range(topological_dimension)]
514
sum_ = format["grouping"](format["add"](inner))
515
value = format["multiply"]([format["inverse"](format["det(J)"]("")), sum_])
516
name = f_component(format["reference derivatives"], _matrix_index(i, f_r, format["num derivatives"]))
517
lines += [f_assign(name, value)]
518
elif mapping == "covariant piola":
519
lines += ["", Indent.indent(format["comment"]\
520
("Using covariant Piola transform to map values back to the physical element"))]
521
# Get temporary values before mapping.
522
lines += [format["const float declaration"](Indent.indent(f_tmp(i)),\
523
f_component(format["reference derivatives"], _matrix_index(i, f_r, format["num derivatives"]))) for i in range(num_components)]
524
# Create names for inner product.
525
topological_dimension = data["topological_dimension"]
526
basis_col = [f_tmp(j) for j in range(topological_dimension)]
527
for i in range(num_components):
528
# Create inverse of Jacobian.
529
inv_jacobian_column = [format["transform"]("JINV", j, i, None) for j in range(topological_dimension)]
531
# Create inner product of basis values and inverse of Jacobian.
532
inner = [format["multiply"]([inv_jacobian_column[j], basis_col[j]]) for j in range(topological_dimension)]
533
value = format["grouping"](format["add"](inner))
534
name = f_component(format["reference derivatives"], _matrix_index(i, f_r, format["num derivatives"]))
535
lines += [f_assign(name, value)]
537
error("Unknown mapping: %s" % mapping)
539
# Generate loop over number of derivatives.
540
code += generate_loop(lines, loop_vars, Indent, format)
544
def _transform_derivatives(data, sum_value_dim, Indent, format):
549
# Prefetch formats to speed up code generation
550
f_loop = format["loop"]
551
f_num_derivatives = format["num derivatives"]
552
f_derivatives = format["reference derivatives"]
553
f_values = format["argument values"]
554
f_mul = format["mul"]
555
f_add = format["add"]
556
f_component = format["component"]
557
f_transform = format["transform matrix"]
559
code += [Indent.indent(format["comment"]("Transform derivatives back to physical element"))]
561
code += [Indent.indent(f_loop("row", 0, f_num_derivatives))]
562
code += [Indent.indent(format["block begin"])]
563
# Increase indentation
565
code += [Indent.indent(f_loop("col", 0, f_num_derivatives))]
566
code += [Indent.indent(format["block begin"])]
567
# Increase indentation
570
# Get number of components, change for tensor valued elements.
571
shape = data["value_shape"]
574
elif len(shape) == 1:
575
num_components = shape[0]
577
error("Tensor valued elements are not supported yet: %s" % data["family"])
579
# Compute offset in array values if any
580
for i in range(num_components):
581
access_name = _matrix_index(sum_value_dim + i, "row", f_num_derivatives)
582
name = f_component(f_values, access_name)
583
# if (sum_value_dim + i == 0):
584
# name = f_component(f_values, "row")
585
# elif (sum_value_dim + i == 1):
586
# name = f_component(f_values, f_add([f_num_derivatives, "row"]))
588
# offset_name = f_mul(["%d" %(sum_value_dim + i), f_num_derivatives])
589
# name = f_component(f_values, f_add([offset_name, "row"]))
590
access_val = _matrix_index(i, "col", f_num_derivatives)
591
value = f_mul([f_component(f_transform, ["row","col"]),\
592
f_component(f_derivatives, access_val)])
594
# value = f_mul([f_component(f_transform, ["row","col"]),\
595
# f_component(f_derivatives, "col")])
597
# value = f_mul([f_component(f_transform, ["row","col"]),\
598
# f_component(f_derivatives, f_add([f_num_derivatives, "col"]))])
600
# offset_value = f_mul(["%d" %i, f_num_derivatives])
602
# value = f_mul([f_component(f_transform, ["row","col"]),\
603
# f_component(f_derivatives, f_add([offset_value, "col"]))])
606
code += [Indent.indent(format["iadd"](name, value))]
608
# Decrease indentation
610
code += [Indent.indent(format["block end"])]
612
# Decrease indentation
615
code += [Indent.indent(format["block end"])]
619
def _delete_pointers(data, Indent, format):
620
"Delete the pointers to arrays."
623
f_r = format["free indices"][0]
625
code += ["", Indent.indent(format["comment"]("Delete pointer to array of derivatives on FIAT element"))]
626
code += [Indent.indent(format["delete pointer"](format["reference derivatives"], "")), ""]
628
code += [Indent.indent(format["comment"]("Delete pointer to array of combinations of derivatives and transform"))]
629
loop_vars = [(f_r, 0, format["num derivatives"])]
630
lines = [format["delete pointer"](format["derivative combinations"], format["component"]("", f_r))]
631
lines += [format["delete pointer"](format["transform matrix"], format["component"]("", f_r))]
632
code += generate_loop(lines, loop_vars, Indent, format)
634
code += [Indent.indent(format["delete pointer"](format["derivative combinations"], ""))]
635
code += [Indent.indent(format["delete pointer"](format["transform matrix"], ""))]