This solves the quadratic scaling of DofMapBuilder::reorder_local;
however there must be something else that is quadratic during the
computation of the matrix sparsity pattern (see example belows with
poisson eq.):
'Real' spaces are a problem and shouldn't be used for large problems.
They lead to dense rows in the matrix. Sparsity pattern construction
and sparse matrix insertion are designed/optimised for sparse rows.
'Real' spaces can be extremely handy in structural mechanics.
They may lead to dense rows in the matrix; however, this is not
necessarily true.
Besides, the bilinear forms below are built without actually using the
lagrange multipliers: the corresponding equations are left
unspecified.
I still have to look at how the sparsity pattern is built in dolfin;
however, I don't see how a lagrange multiplier could possibly
lead to a quadratic behaviour in a "tradizional" code that computes the
sparsity pattern by looping over all the elements.
Marco
Garth
#========================================
from dolfin import *
import time
import numpy as np
import matplotlib.pyplot as plt
set_log_level(DEBUG)
parameters['reorder_dofs_serial'] = False
for n in range(25, 33):
mesh = UnitCubeMesh(n, n, n)
#mesh = UnitSquareMesh(10*n, 10*n)
P1 = FunctionSpace(mesh, 'CG', 1)
R = FunctionSpace(mesh, 'R', 0)
# Measure time for creation of P1*R space
t = time.time()
V = P1*R
t = time.time() - t
PP = P1*P1
(uu, ul) = TrialFunctions(V)
(vu, vl) = TestFunctions(V)
(pu, pp) = TrialFunctions(PP)
(pv, pq) = TestFunctions(PP)
u = TrialFunction(P1)
v = TestFunction(P1)
am = inner(grad(uu), grad(vu))*dx
a = inner(grad(pu), grad(pv))*dx
af = inner(grad(u), grad(v))*dx
#mixed cg-real
t1 = time.time()
AM = assemble(am)
t1 = time.time() - t1
#mixed cg-cg
t2 = time.time()
A = assemble(a)
t2 = time.time() - t2
#cg
t3 = time.time()
AF = assemble(af)
t3 = time.time() - t3
# Report
print 'DOFs:', V.dim(),\
'T/D:', t/V.dim(),\
'T1/D:', t1/P1.dim(),\
'T2/D:', t2/P1.dim(),\
'T3/D:', t3/P1.dim()
plt.loglog(V.dim(), t1, 'o')
plt.loglog(V.dim(), t2, 'x')
plt.loglog(V.dim(), t3, '*')
#
# Plot some quadratic function for comparison
xl = plt.xlim()
yl = plt.ylim()
x = np.linspace(*xl)
y = map(lambda x:yl[0]/(xl[0]*xl[0])*x*x, x)
plt.loglog(x, y)
#
# # Show plot
plt.show()
#=============================================
Marco
.
--
Dr. Marco Morandini, Ph.D.
Politecnico di Milano
Dip. di Scienze e Tecnologie Aerospaziali
Via La Masa, 34 | phone: (+39) 02 2399 8362
20156 MILANO - Italy | fax: (+39) 02 2399 8334
-----------------------------------------------------------
Dr. Marco Morandini, Ph.D.
Dipartimento di Scienze e Tecnologie Aerospaziali
Politecnico di Milano
tel: (+39) 02 2399 8362
fax: (+39) 02 2399 8334
-----------------------------------------------------------
Dr. Marco Morandini, Ph.D. | mail: [email protected]
Via J. F. Kennedy 20 | cell: (+39) 340 7956882
20097 S. Donato M.se (MI) - Italy | home: (+39) 02 54100118
-----------------------------------------------------------
Dr. Marco Morandini, Ph.D.
via J. F. Kennedy 20
20097 S. Donato M.se (MI) - Italy
mail: [email protected]
cell: (+39) 340 7956882
home: (+39) 02 54100118
-----------------------------------------------------------
Please avoid sending me Word or PowerPoint attachments.
See http://www.gnu.org/philosophy/no-word-attachments.html
_______________________________________________
fenics mailing list
[email protected]
http://fenicsproject.org/mailman/listinfo/fenics
