On 6 June 2011 12:33, Anders Logg <l...@simula.no> wrote: > On Mon, Jun 06, 2011 at 10:17:27AM -0000, Garth Wells wrote: >> On 06/06/11 11:05, Martin Sandve Alnæs wrote: >> > On 6 June 2011 11:54, Garth Wells <785...@bugs.launchpad.net> wrote: >> >> On 06/06/11 10:41, Martin Sandve Alnæs wrote: >> >>> On 31 May 2011 00:24, Anders Logg <l...@simula.no> wrote: >> >>>> On Mon, May 30, 2011 at 09:53:42PM -0000, Martin Sandve Alnæs wrote: >> >>>>> There's two, don't remember what they do: >> >>>>> def estimate_max_polynomial_degree(e, default_degree=1): >> >>>>> def estimate_total_polynomial_degree(e, default_degree=1): >> >>>>> in algorithms/transformations.py (should rather be in analysis.py I >> >>>>> guess). >> >>>>> >> >>>>> ** Changed in: dolfin >> >>>>> Status: New => Invalid >> >>>> >> >>>> And those include spatial coordinates? >> >>> >> >>> Turns out they didn't. Just checked the code. >> >>> But it was easy to add. I'm commiting changes >> >>> to estimate_total_polynomial_degree now which >> >>> incorporate the spatial degree. Maybe this should >> >>> be used for assembling rhs and functionals, while >> >>> looking at elements are enough for the bilinear form? >> >>> >> >>> PyDOLFIN could do something like >> >>> >> >>> d = estimate_total_polynomial_degree(expr) >> >>> d = max(d, 1) >> >>> d = min(d, 8) >> >>> >> >>> to limit the degree to some reasonable range in cases such as >> >>> expr = sin(x**5)*cos(y**5) >> >>> which would lead to a degree of (5+2)+(5+2)=14 with the current >> >>> heuristics. >> >>> Look at the code and tests in the last commit for more details, it's >> >>> quite short. >> >>> >> >> >> >> We have the same issue of order blow-out for problems with lots of >> >> coefficients. I'm therefore inclined not to include any heuristics, and >> >> leave it up to the user. >> > >> > Do you mean we should actually crash and burn with this line? >> > f = assemble(sin(triangle.x[0]), mesh=mesh) >> > With the current heuristic this will give degree 3, >> > 1 from x and +2 from sin. >> > >> >> We obviously need an approach for functions from non-polynomial spaces. >> What I'm not inclined towards is arbitrary thresholds for integrating >> polynomial products. > > What if we by default set that threshold to the maximal element degree > + k, where k is say 1? It would be enough to retain expected > convergence. > > -- > Anders
That's definitely not enough (e.g. L = x**8*v*dx) and would be basically the same as throwing away the entire total polynomial degree algorithm. It might be necessary/desirable to use different rules for bilinear forms and linear forms/functionals? Martin -- You received this bug notification because you are a member of DOLFIN Team, which is subscribed to DOLFIN. https://bugs.launchpad.net/bugs/785874 Title: Projection of x is not accurate Status in DOLFIN: Invalid Bug description: I've tested that projecting x works without the scaling bug that was just fixed, using dimensions 1,2,3 and both DG and CG from 0 to 3 degrees. I print the max and min values of the vector of the projection function, and the values are _close_ to 0 and 1 but not to machine precision. The script is below. There's up to 2.7% error in the 3D case. Is the projection form integrated accurately enough? All but the DG0 function space should be capable of representing x exactly. Not sure if this is a dolfin or ffc bug. from dolfin import * def mcx(dim): if dim == 1: mesh = UnitInterval(20) cell = interval x = cell.x if dim == 2: mesh = UnitSquare(20, 20) cell = triangle x = cell.x[0] if dim == 3: mesh = UnitCube(20, 20, 20) cell = tetrahedron x = cell.x[0] return mesh, cell, x for dim in range(1, 4): mesh, cell, x = mcx(dim) minval, maxval = 1.0, 0.0 #print dim, "DG" for degree in range(3): V = FunctionSpace(mesh, "DG", degree) u = project(x, V) #print dim, degree, u.vector().min(), u.vector().max() minval = min(minval, u.vector().min()) maxval = max(maxval, u.vector().max()) #print dim, "CG" for degree in range(1, 4): V = FunctionSpace(mesh, "CG", degree) u = project(x, V) #print dim, degree, u.vector().min(), u.vector().max() minval = min(minval, u.vector().min()) maxval = max(maxval, u.vector().max()) print minval, maxval _______________________________________________ Mailing list: https://launchpad.net/~dolfin Post to : dolfin@lists.launchpad.net Unsubscribe : https://launchpad.net/~dolfin More help : https://help.launchpad.net/ListHelp
