Screenshot instructions:
Windows
Mac
Red Hat Linux
Ubuntu
Click URL instructions:
Rightclick on ad, choose "Copy Link", then paste here →
(This may not be possible with some types of ads)
You can subscribe to this list here.
2010 
_{Jan}

_{Feb}

_{Mar}

_{Apr}

_{May}

_{Jun}

_{Jul}

_{Aug}

_{Sep}

_{Oct}
(1) 
_{Nov}
(1) 
_{Dec}


2011 
_{Jan}

_{Feb}

_{Mar}
(1) 
_{Apr}

_{May}
(1) 
_{Jun}
(9) 
_{Jul}
(2) 
_{Aug}
(6) 
_{Sep}
(11) 
_{Oct}
(15) 
_{Nov}
(4) 
_{Dec}
(9) 
2012 
_{Jan}
(7) 
_{Feb}
(14) 
_{Mar}
(11) 
_{Apr}

_{May}

_{Jun}
(3) 
_{Jul}

_{Aug}

_{Sep}

_{Oct}

_{Nov}
(3) 
_{Dec}

2013 
_{Jan}
(8) 
_{Feb}

_{Mar}
(1) 
_{Apr}
(3) 
_{May}

_{Jun}

_{Jul}
(3) 
_{Aug}
(2) 
_{Sep}

_{Oct}

_{Nov}

_{Dec}
(4) 
2014 
_{Jan}

_{Feb}

_{Mar}

_{Apr}

_{May}
(2) 
_{Jun}
(2) 
_{Jul}

_{Aug}

_{Sep}
(5) 
_{Oct}
(3) 
_{Nov}
(1) 
_{Dec}
(4) 
2015 
_{Jan}
(5) 
_{Feb}
(10) 
_{Mar}
(12) 
_{Apr}
(14) 
_{May}
(26) 
_{Jun}
(7) 
_{Jul}

_{Aug}

_{Sep}
(4) 
_{Oct}

_{Nov}

_{Dec}

2016 
_{Jan}
(1) 
_{Feb}
(9) 
_{Mar}
(3) 
_{Apr}
(3) 
_{May}
(12) 
_{Jun}
(2) 
_{Jul}
(16) 
_{Aug}
(4) 
_{Sep}
(18) 
_{Oct}
(3) 
_{Nov}
(4) 
_{Dec}
(10) 
2017 
_{Jan}

_{Feb}

_{Mar}
(9) 
_{Apr}
(14) 
_{May}
(12) 
_{Jun}
(17) 
_{Jul}

_{Aug}
(3) 
_{Sep}
(1) 
_{Oct}
(1) 
_{Nov}

_{Dec}

2018 
_{Jan}
(7) 
_{Feb}
(3) 
_{Mar}
(2) 
_{Apr}
(2) 
_{May}

_{Jun}

_{Jul}

_{Aug}

_{Sep}

_{Oct}

_{Nov}

_{Dec}

From: Rafael Vazquez <vazquez@im...>  20180419 08:00:52

Dear Tony, to work with orthotropic materials, you will first need to modify the operator functions, like op_su_ev, in such a way that the physical parameters are not scalars, but tensors. You will also need to modify the lambda_lame and mu_lame functions, in such a way that the output has the correct size. I don't have any experience with orthotropic materials, so I don't know if there are any other necessary modifications. Best, Rafa On 18. 04. 18 16:37, TonyZ wrote: > Dear Rafa > > I am working on the simulation of orthotropic material. For the 2D > linear elasticity problem, is it possible to use GeoPDEs > to solve orthotropic problems? how can I modify following code? > > E = [25.5e9 25.5e9 28e9 27e9 28e9 ]; > nu = [0.167 0.167 0.200 0.200 0.200]; > problem_data.lambda_lame = @(x, y, iptc) > ((nu(iptc)*E(iptc))/((1+nu(iptc))*(12*nu(iptc))) * ones (size (x))); > problem_data.mu_lame = @(x, y, iptc) (E(iptc)/(2*(1+nu(iptc))) * ones > (size (x))); > > Any suggestions are gratefully accepted. > > Yours Sincerely, > Tony > > 
From: TonyZ <tju_zmx@12...>  20180418 14:37:24

Dear Rafa I am working on the simulation of orthotropic material. For the 2D linear elasticity problem, is it possible to use GeoPDEs to solve orthotropic problems? how can I modify following code? E = [25.5e9 25.5e9 28e9 27e9 28e9 ]; nu = [0.167 0.167 0.200 0.200 0.200]; problem_data.lambda_lame = @(x, y, iptc) ((nu(iptc)*E(iptc))/((1+nu(iptc))*(12*nu(iptc))) * ones (size (x))); problem_data.mu_lame = @(x, y, iptc) (E(iptc)/(2*(1+nu(iptc))) * ones (size (x))); Any suggestions are gratefully accepted. Yours Sincerely, Tony 
From: Rafael Vazquez <vazquez@im...>  20180320 13:04:17

Dear Ganesh, the geopdes_tsplines package was intended to work with the files from the Tspline plugin for Rhino, and it was tested with Rhino 4.0. There was no further development of the package, and now it is not part of the latest release of GeoPDEs. Apart from Rhino, I am not aware of any other code generating an analogous file, but I have not worked with Tsplines for more than five years now. Best regards, Rafael On 20. 03. 18 13:31, Ganesh Diwan wrote: > Dear Developers, > > I am trying my hands on geopdes_tsplines function available in the library. > It appears that one always has to use the Tspline files generated by > Autodesk's plugin. I have been trying to get Tspline plugin (version 4 ) > work with Rhino version 5 (demo version) without any success so far. Is > there a possibility one can generate the tspline files for canonical > geometries (in the format the function read_bezier_extraction needs them ) > without using Rhino? > > Much appreciate any help / comments. Thanks. > > Ganesh Diwan. >  > Check out the vibrant tech community on one of the world's most > engaging tech sites, Slashdot.org! http://sdm.link/slashdot > _______________________________________________ > Geopdesusers mailing list > Geopdesusers@... > https://lists.sourceforge.net/lists/listinfo/geopdesusers 
From: Ganesh Diwan <gcdiwan@gm...>  20180320 12:31:57

Dear Developers, I am trying my hands on geopdes_tsplines function available in the library. It appears that one always has to use the Tspline files generated by Autodesk's plugin. I have been trying to get Tspline plugin (version 4 ) work with Rhino version 5 (demo version) without any success so far. Is there a possibility one can generate the tspline files for canonical geometries (in the format the function read_bezier_extraction needs them ) without using Rhino? Much appreciate any help / comments. Thanks. Ganesh Diwan. 
From: Wilfred van Rooijen <wvanrooijen@ya...>  20180227 12:30:14

Well....... About the geometry: if you have no experience with NURBS geometry, you must get a copy of "The NURBS Book", by Tiller & Piegl. It is likely that your university library has a copy, probably in the engineering department, or architecture, any place where CAD is used. GeoPDEs has many functions to create simple NURBS geometries. Start Octave, then > pkg load geopdes> help nrb+TAB+TAB where +TAB+TAB means "hit the TAB key twice". All the functions starting with "nrb" will be displayed, including things like nrb4surf (a 2D rectangular domain), nrbcirc (part of a circle), nrbcylinder, etc. The helpfunction will then show you how to proceed. About the problem to solve, start out by writing down the equations you want to solve, then determine the necessary operators for the Galerkin method, and use the operators from GeoPDEs, like op_u_v_tp, op_vel_gradu_tp, op_gradv_gradu_tp, etc. I programmed a solver for the neutron diffusion equation in GeoPDEs last year. I think the heat diffusion (conduction) is very similar, although I only work on static problems (no time dependence). Wilfred On Tuesday, February 27, 2018 6:20 PM, Rafael Vazquez <vazquez@...> wrote: Dear Jan, regarding the simulation without a proper geometry. The simplest thing you can do is to define the geometry as the identity. You can either do it with the NURBS toolbox, or just define the functions for the parametrization and the derivative, to save computational time (see the example ex_article_section_512.m, and the functions ring_polar_map and ring_polar_map_der). This will also allow you to generalize to the case of having a geometry. The alternative, to completely avoid the definition of the geometry, would be to change the definition of the msh and space objects, and change also the operators inside the class (and some other functions) to compute the basis functions in the parametric domain with sp_evaluate_col_param. Trust me, you don't want to do that. ;) About the time scheme, let's see if other users can share their code with you. I have an old implicit Euler scheme implementation (it should be somewhere in my laptop), but I don't know how useful it will be if you want to implement a Galerkin method in time. Best regards, Rafa On 27. 02. 18 09:13, Jan Pospíšil wrote: > Hello all from the GeoPDEs community. > > At the IGA conference in Pavia last year I made a note about GeoPDEs that I finally got some time to install and now I would like to test if it is suitable for modification in order to solve my type of problems (partial integrodifferential equations  PIDEs). Rafael Vazques directed me to this mailing list and I would like to ask you for your help. I think that I already got used to the syntax of the functions in the NURBS toolbox GeoPDEs is using and now I have started with the first examples. > > I would like to start with a problem, where I do not have a proper NURBS geometry, say with a simple example with finite rectangular domain in one or two spatial and one temporal dimension that I would like to solve using FEM with NURBS (or at least Bspline) basis. Say for example the (parabolic) heat equation in one spatial and one temporal dimension, where I know the exact solution that I could use for example for comparison and to do for example some optimisation in the number of control points and in the basis order in order to get a for example a prescribed l2 error etc. > > Has somebody out there already tried GeoPDEs for solving parabolic PDEs? If I am not mistaken, all provided examples are for elliptic equations only, am I right? A sufficient working example could use FEM only in spatial dimensions and for example a fully implicit time stepping scheme. > > Several your examples show how to use IGA for problems on unit square domains. Do I understand right that a geometry is always required in GeoPDEs in order to have a suitable “mesh”? May I ask you how would you solve the above mentioned problem where there is no proper geometry? A fully working example would be highly appreciated. > > Kind regards > Jan > >  > Jan Pospisil, Ph.D. email: honik@... > University of West Bohemia phone: (+420) 377632675 > Department of Mathematics fax: (+420) 377632602 > Plzen, Czech Republic address: Univerzitni 8, 306 14 > > >  > Check out the vibrant tech community on one of the world's most > engaging tech sites, Slashdot.org! http://sdm.link/slashdot > _______________________________________________ > Geopdesusers mailing list > Geopdesusers@... > https://lists.sourceforge.net/lists/listinfo/geopdesusers  Check out the vibrant tech community on one of the world's most engaging tech sites, Slashdot.org! http://sdm.link/slashdot _______________________________________________ Geopdesusers mailing list Geopdesusers@... https://lists.sourceforge.net/lists/listinfo/geopdesusers 
From: Rafael Vazquez <vazquez@im...>  20180227 09:20:45

Dear Jan, regarding the simulation without a proper geometry. The simplest thing you can do is to define the geometry as the identity. You can either do it with the NURBS toolbox, or just define the functions for the parametrization and the derivative, to save computational time (see the example ex_article_section_512.m, and the functions ring_polar_map and ring_polar_map_der). This will also allow you to generalize to the case of having a geometry. The alternative, to completely avoid the definition of the geometry, would be to change the definition of the msh and space objects, and change also the operators inside the class (and some other functions) to compute the basis functions in the parametric domain with sp_evaluate_col_param. Trust me, you don't want to do that. ;) About the time scheme, let's see if other users can share their code with you. I have an old implicit Euler scheme implementation (it should be somewhere in my laptop), but I don't know how useful it will be if you want to implement a Galerkin method in time. Best regards, Rafa On 27. 02. 18 09:13, Jan Pospíšil wrote: > Hello all from the GeoPDEs community. > > At the IGA conference in Pavia last year I made a note about GeoPDEs that I finally got some time to install and now I would like to test if it is suitable for modification in order to solve my type of problems (partial integrodifferential equations  PIDEs). Rafael Vazques directed me to this mailing list and I would like to ask you for your help. I think that I already got used to the syntax of the functions in the NURBS toolbox GeoPDEs is using and now I have started with the first examples. > > I would like to start with a problem, where I do not have a proper NURBS geometry, say with a simple example with finite rectangular domain in one or two spatial and one temporal dimension that I would like to solve using FEM with NURBS (or at least Bspline) basis. Say for example the (parabolic) heat equation in one spatial and one temporal dimension, where I know the exact solution that I could use for example for comparison and to do for example some optimisation in the number of control points and in the basis order in order to get a for example a prescribed l2 error etc. > > Has somebody out there already tried GeoPDEs for solving parabolic PDEs? If I am not mistaken, all provided examples are for elliptic equations only, am I right? A sufficient working example could use FEM only in spatial dimensions and for example a fully implicit time stepping scheme. > > Several your examples show how to use IGA for problems on unit square domains. Do I understand right that a geometry is always required in GeoPDEs in order to have a suitable “mesh”? May I ask you how would you solve the above mentioned problem where there is no proper geometry? A fully working example would be highly appreciated. > > Kind regards > Jan > >  > Jan Pospisil, Ph.D. email: honik@... > University of West Bohemia phone: (+420) 377632675 > Department of Mathematics fax: (+420) 377632602 > Plzen, Czech Republic address: Univerzitni 8, 306 14 > > >  > Check out the vibrant tech community on one of the world's most > engaging tech sites, Slashdot.org! http://sdm.link/slashdot > _______________________________________________ > Geopdesusers mailing list > Geopdesusers@... > https://lists.sourceforge.net/lists/listinfo/geopdesusers 
From: Jan Pospíšil <honik@km...>  20180227 08:13:33

Hello all from the GeoPDEs community. At the IGA conference in Pavia last year I made a note about GeoPDEs that I finally got some time to install and now I would like to test if it is suitable for modification in order to solve my type of problems (partial integrodifferential equations  PIDEs). Rafael Vazques directed me to this mailing list and I would like to ask you for your help. I think that I already got used to the syntax of the functions in the NURBS toolbox GeoPDEs is using and now I have started with the first examples. I would like to start with a problem, where I do not have a proper NURBS geometry, say with a simple example with finite rectangular domain in one or two spatial and one temporal dimension that I would like to solve using FEM with NURBS (or at least Bspline) basis. Say for example the (parabolic) heat equation in one spatial and one temporal dimension, where I know the exact solution that I could use for example for comparison and to do for example some optimisation in the number of control points and in the basis order in order to get a for example a prescribed l2 error etc. Has somebody out there already tried GeoPDEs for solving parabolic PDEs? If I am not mistaken, all provided examples are for elliptic equations only, am I right? A sufficient working example could use FEM only in spatial dimensions and for example a fully implicit time stepping scheme. Several your examples show how to use IGA for problems on unit square domains. Do I understand right that a geometry is always required in GeoPDEs in order to have a suitable “mesh”? May I ask you how would you solve the above mentioned problem where there is no proper geometry? A fully working example would be highly appreciated. Kind regards Jan  Jan Pospisil, Ph.D. email: honik@... University of West Bohemia phone: (+420) 377632675 Department of Mathematics fax: (+420) 377632602 Plzen, Czech Republic address: Univerzitni 8, 306 14 
From: Rafael Vazquez <vazquez@im...>  20180126 19:07:11

Hello Wilfred, if I understood correctly, the algorithms that you are looking for are implemented in the function "nrbeval", to calculate a point, and in the functions "nrbderiv" and "nrbdeval" for the derivatives. They should be already implemented for volumes. Please, check if this is correct. By the way, in GeoPDEs you can directly use the function handles stored in the "geometry" structure: geometry.map and geometry.map_der (for convenience they are also stored in msh). And if it is not what you are looking for, please try to detail your question, and we will find a solution. :) Best, Rafa On 26. 01. 18 12:39, Wilfred van Rooijen wrote: > Hello Rafael, > > The NURBS book contains algorithms for 2D surface, to calculate a > "point" (i.e. determine x(u0, v0) and y(u0, v0)) as well as the > derivatives of the surface at that point. For 3D NURBS volumes, > GeoPDEs must contain the same algorithms, but generalized to 3D. > > OK, I will reread Piegl & Tiller and try to develop the necessary > algorithms. I will then use GeoPDEs to check my algorithms. > > Cheers, > Wilfred van Rooijen > > > On Wednesday, January 24, 2018 6:32 AM, Rafael Vazquez > <vazquez@...> wrote: > > > Hello Wilfred, > sorry, I'm in travel and I don't have access to the NURBS book. Can > you please detail which functions of GeoPDEs (or the NURBS package) > correspond to those algorithms (if any), and the algorithm number in > the NURBS book? > > In any case, the volume functions that I already implemented for the > NURBS package were the ones that were easy to generalize from surfaces > to volumes. If the the functions that you need are not implemented for > volumes, probably you will need to implement them (I can try to give > you some help). As far as I know, there is no reference for the > algorithms of volumes/trivariates, but I hope that calculating a point > and the derivatives can be "easily" generalized. > > Best, > Rafa > > > On 22. 01. 18 06:49, Wilfred van Rooijen wrote: > Hello Rafael, others, > > My student is currently working on the development of an IGAbased > code for neutron transport theory. We have a 2D version and I am now > investigation a 3D version. > > For the 2D version, I have used several algorithms from Piegl & Tiller > "The NURBS Book": > >  algorithm to calculate a "point" on a surface (SURFACEPOINT) >  algorithm to calculate a point on a surface and the derivatives > (RATSURFDERIVS) > > The NURBS book does not have the 3D algorithms. Can you recommend a > good source for the 3D algorithms, or should I go and read the NURBS > coding in GeoPDEs? > > Thanks, > Wilfred > > > > 
From: Wilfred van Rooijen <wvanrooijen@ya...>  20180126 11:39:54

Hello Rafael, The NURBS book contains algorithms for 2D surface, to calculate a "point" (i.e. determine x(u0, v0) and y(u0, v0)) as well as the derivatives of the surface at that point. For 3D NURBS volumes, GeoPDEs must contain the same algorithms, but generalized to 3D. OK, I will reread Piegl & Tiller and try to develop the necessary algorithms. I will then use GeoPDEs to check my algorithms. Cheers,Wilfred van Rooijen On Wednesday, January 24, 2018 6:32 AM, Rafael Vazquez <vazquez@...> wrote: Hello Wilfred, sorry, I'm in travel and I don't have access to the NURBS book. Can you please detail which functions of GeoPDEs (or the NURBS package) correspond to those algorithms (if any), and the algorithm number in the NURBS book? In any case, the volume functions that I already implemented for the NURBS package were the ones that were easy to generalize from surfaces to volumes. If the the functions that you need are not implemented for volumes, probably you will need to implement them (I can try to give you some help). As far as I know, there is no reference for the algorithms of volumes/trivariates, but I hope that calculating a point and the derivatives can be "easily" generalized. Best, Rafa On 22. 01. 18 06:49, Wilfred van Rooijen wrote: Hello Rafael, others, My student is currently working on the development of an IGAbased code for neutron transport theory. We have a 2D version and I am now investigation a 3D version. For the 2D version, I have used several algorithms from Piegl & Tiller "The NURBS Book":  algorithm to calculate a "point" on a surface (SURFACEPOINT)  algorithm to calculate a point on a surface and the derivatives (RATSURFDERIVS) The NURBS book does not have the 3D algorithms. Can you recommend a good source for the 3D algorithms, or should I go and read the NURBS coding in GeoPDEs? Thanks, Wilfred 
From: Rafael Vazquez <vazquez@im...>  20180123 21:56:47

Hello Wilfred, sorry, I'm in travel and I don't have access to the NURBS book. Can you please detail which functions of GeoPDEs (or the NURBS package) correspond to those algorithms (if any), and the algorithm number in the NURBS book? In any case, the volume functions that I already implemented for the NURBS package were the ones that were easy to generalize from surfaces to volumes. If the the functions that you need are not implemented for volumes, probably you will need to implement them (I can try to give you some help). As far as I know, there is no reference for the algorithms of volumes/trivariates, but I hope that calculating a point and the derivatives can be "easily" generalized. Best, Rafa On 22. 01. 18 06:49, Wilfred van Rooijen wrote: > Hello Rafael, others, > > My student is currently working on the development of an IGAbased > code for neutron transport theory. We have a 2D version and I am now > investigation a 3D version. > > For the 2D version, I have used several algorithms from Piegl & Tiller > "The NURBS Book": > >  algorithm to calculate a "point" on a surface (SURFACEPOINT) >  algorithm to calculate a point on a surface and the derivatives > (RATSURFDERIVS) > > The NURBS book does not have the 3D algorithms. Can you recommend a > good source for the 3D algorithms, or should I go and read the NURBS > coding in GeoPDEs? > > Thanks, > Wilfred > 
From: Wilfred van Rooijen <wvanrooijen@ya...>  20180122 05:50:15

Hello Rafael, others, My student is currently working on the development of an IGAbased code for neutron transport theory. We have a 2D version and I am now investigation a 3D version. For the 2D version, I have used several algorithms from Piegl & Tiller "The NURBS Book":  algorithm to calculate a "point" on a surface (SURFACEPOINT) algorithm to calculate a point on a surface and the derivatives (RATSURFDERIVS) The NURBS book does not have the 3D algorithms. Can you recommend a good source for the 3D algorithms, or should I go and read the NURBS coding in GeoPDEs? Thanks,Wilfred 
From: Sander Dedoncker <sander.dedoncker@ku...>  20180108 14:35:08

Dear Rafael, I'm sorry, there was indeed a mistake on my part. I wasn't careful enough in cleaning my Matlab path, and didn't notice that an older version of nrbextract was being called. Setting this right fixed the error for me as well. Thanks a lot for the quick reply and my apologies for the false alarm! Best regards, Sander Oorspronkelijk bericht Van: Rafael Vazquez [mailto:vazquez@...] Verzonden: maandag 8 januari 2018 14:40 Aan: Sander Dedoncker <sander.dedoncker@...>; geopdesusers@... Onderwerp: Re: [Geopdesusers] Bug in mp_geo_load? Dear Sander, I was unable to reproduce the error. The code you sent works for me, both in Octave and Matlab, using GeoPDEs version 3.1.0 and the NURBS package version 1.3.13. Are you sure to be using the correct version of the NURBS package? The function nrbextract was generalized for curves in version 1.3.11, and the error message you get is no longer present. Regards, Rafael 
From: Rafael Vazquez <vazquez@im...>  20180108 14:04:56

Dear Sander, I was unable to reproduce the error. The code you sent works for me, both in Octave and Matlab, using GeoPDEs version 3.1.0 and the NURBS package version 1.3.13. Are you sure to be using the correct version of the NURBS package? The function nrbextract was generalized for curves in version 1.3.11, and the error message you get is no longer present. Regards, Rafael 
From: Sander Dedoncker <sander.dedoncker@ku...>  20180104 15:37:23

Dear, I am trying to use the routines in GeoPDEs to determine the interfaces and boundaries in a multipatch NURBS geometry. However I seem to encounter a bug when calling 'mp_geo_load'. Consider the following simple example: generate an array of two NURBS structures representing adjacent squares, using cpts1 = cat(3,[0 1; 0 0; 0 0; 1 1],[0 1; 1 1; 0 0; 1 1]); cpts2 = cpts1; cpts2(1,:) = cpts2(1,:)+1; knots = {[0 0 1 1],[0 0 1 1]}; patch1 = nrbmak(cpts1,knots); patch2 = nrbmak(cpts2,knots); mpgeom = [patch1;patch2]; When I now call mp_geo_load, [geometry, boundaries, interfaces, subdomains, boundary_interfaces] = mp_geo_load (mpgeom); the code errors as follows: Error using nrbextract (line 47) The boundary information is only extracted for NURBS surfaces or volumes Error in nrbmultipatch (line 79) nrb_faces1 = nrbextract (nurbs(i1)); Error in mp_geo_load (line 192) boundary_interfaces = nrbmultipatch (bnd_nurbs); The problem occurs in GeoPDEs version 3.1.0 with NURBS Toolbox version 1.3.13. From the documentation of mp_geo_load, I understand that this case should be easily treated. Is that correct? Thanks and best regards, Sander Dedoncker 
From: Rafael Vazquez <vazquez@im...>  20171012 12:57:43

Dear GeoPDEs users, we have included in the NURBS package two new functions (nrbspheretiling and nrbspheretile) to generate different multipatch parametrizations of the sphere. The functions are available in the NURBS repository in Octaveforge, and will be released with the new version of the package. Thanks to Sander Dedoncker from KU Leuven, and his collaborators, for providing them. Best regards, Rafa 
From: Rafael Vazquez <vazquez@im...>  20170901 12:59:24

Dear Léo, the corners are not stored by default, but you should be able to compute them quite easily. You only need to create an auxiliary mesh object, with the points where you want to compute, and the breakpoints between elements that you want (can be the knots, or whatever you want). It is the same that is done in msh_side_from_interior, but for the corners you will need to do it for all the directions. After generating the auxiliary mesh, you can call the constructor of the space that you already had space_corners = space.constructor (msh_auxiliary); When you compute using this new space object, you will obtain the values of the basis functions at the corner points, instead of the quadrature points. Regarding the axisymmetric problem, yes, I modified all the operators you mention. I can only suggest to try the single patch problem with known analytical solution. That should help you to understand where the error is coming from. Best regards, Rafa On 31. 08. 17 18:30, Friedrich, L.A.J. wrote: > Dear Rafa, > > Thank you for your answers, it solved (almost completely) my problems for Q1 and Q2. > > I wanted to ask you about corner points in GeoPDE. > I understood that the natural quadrature nodes are only in the interior of the domains, and that we can build a mesh on the interfaces between patches, to obtain quadrature nodes on that boundary, that allow to impose boundary conditions or calculate integrals. But I am under the impression that corners points are absent (in the integral calculation and boundary imposition) therefore leading to small errors/discrepancy. Could you please comment on this point, and tell me if/how it is possible to access those corner points ? I thought the answer would lie in the msh_side.breaks, but I do not know how to construct a msh struct from a cell, and the msh_side.bounday is null []. > > For Q3. > I actually have very large discrepancy (4.5e3%) my solution is by order of magnitude larger that it should. My dimensions are in millimeter, but it looks like changing both the stiffness operators and rhs kind of cancel out the multiplication by r in the integrands, even though it should be done like that (30% discrepancy were actually for a strange test where I computed the stiffness matrix in Carthesian and the rhs in axisymmetric). > Therefore I do not really know where to investigate further, except as you suggested on a single patch domain, as it should have been quite straightforward modification. > To clarify you previous comment, you did modify all the following operators: op_gradu_gradv_tp, op_f_v_tp, and in sp_proj_l2_rchlt : op_u_v_tp right ? not just op_gradu_gradv_tp. > > > Best regards, > > Léo > > >  > Check out the vibrant tech community on one of the world's most > engaging tech sites, Slashdot.org! http://sdm.link/slashdot > _______________________________________________ > Geopdesusers mailing list > Geopdesusers@... > https://lists.sourceforge.net/lists/listinfo/geopdesusers 
From: Friedrich, L.A.J. <L.Friedrich@tu...>  20170831 16:46:10

Dear Rafa, Thank you for your answers, it solved (almost completely) my problems for Q1 and Q2. I wanted to ask you about corner points in GeoPDE. I understood that the natural quadrature nodes are only in the interior of the domains, and that we can build a mesh on the interfaces between patches, to obtain quadrature nodes on that boundary, that allow to impose boundary conditions or calculate integrals. But I am under the impression that corners points are absent (in the integral calculation and boundary imposition) therefore leading to small errors/discrepancy. Could you please comment on this point, and tell me if/how it is possible to access those corner points ? I thought the answer would lie in the msh_side.breaks, but I do not know how to construct a msh struct from a cell, and the msh_side.bounday is null []. For Q3. I actually have very large discrepancy (4.5e3%) my solution is by order of magnitude larger that it should. My dimensions are in millimeter, but it looks like changing both the stiffness operators and rhs kind of cancel out the multiplication by r in the integrands, even though it should be done like that (30% discrepancy were actually for a strange test where I computed the stiffness matrix in Carthesian and the rhs in axisymmetric). Therefore I do not really know where to investigate further, except as you suggested on a single patch domain, as it should have been quite straightforward modification. To clarify you previous comment, you did modify all the following operators: op_gradu_gradv_tp, op_f_v_tp, and in sp_proj_l2_rchlt : op_u_v_tp right ? not just op_gradu_gradv_tp. Best regards, Léo 
From: Rafael Vazquez <vazquez@im...>  20170827 11:59:07

Dear Friedrich, I'm glad to see that you are making good progress with GeoPDEs. I reply to your questions below. On 25. 08. 17 15:04, Friedrich, L.A.J. wrote: > Q1. > > In one simulation, I want to impose a source term on an interface between 2 patches, more precisely on top and bottom boundary of one patch. > Therefore I used an adapted version of op_f_v_tp where I start to compute the interior rhs by column, the interior points are set to 0. Then I wrote > ... > for iside = 3 > msh_side = msh_eval_boundary_side (msh, iside); > sp_side = sp_eval_boundary_side (space, msh_side); > for idim = 1:msh.rdim > x{idim} = reshape (msh_side.geo_map(idim,:,:), ...msh_side.nqn, msh_side.nel); > > end > coeff_top =f.t(x{:})%inherited from additional input funtions > rhs = rhs + op_f_v (sp_side, msh_side, coeff_top); > end > ... > > > > same lines for iside=4, but with coeff_bot as the function. > > > > The thing is in one patch I have 104dofs (nsub=[9,4] and nquad=[4,4], degree =[4,4]) which is the same when I compute the rhs by column, but then, on one boundary I have only 13dofs. > > Therefore how to assemble rhs as they have different dimensions ? The boundary entities contain the information about basis functions that do not vanish on the boundary. You can recover their global/volumetric numbering using space.boundary(side).dofs, exactly as for Neumann boundary conditions (see Section 3.7 of the paper). In the case of a multipatch domain, you will also need to apply space.gnum, to pass from the volumetric number of the patch, to the global multipatch numbering. > > Additional question here : we have quadrature points in the interior of the domain, but also on the boundaries right ? is there a way to plot them on top of the geometry? You can plot the geometry with nrbkntplot(geometry.nurbs). In your example, the boundary quadrature points are stored in sp_side.geo_map. You can plot them like this: X = sp_side.geo_map(1,:,:); Y = sp_side.geo_map(2,:,:); Z = sp_side.geo_map(3,:,:); plot3 (X(:), Y(:), Z(:), 'x') Do not forget to use "hold on", to plot both together. > > > > > Q2. > > In another simulation, I want in a postprocessing part, integrate a quantity (deriving from the gradient of my solution) around the boundary of one (group of) patch. Here is the part of code I start to write for one boundary of one patch > function [FX] = sp_eval_FX (u,space,msh,geometry,pts,mat) > FX=0; > for iptc=1:space.npatch > if ismember([iptc],mat) > for iside = 3 > msh_side = msh_eval_boundary_side (msh.msh_patch{iptc}, iside); > sp_side = sp_eval_boundary_side (space.sp_patch{iptc}, msh_side); > for idim = 1:msh.rdim > x{idim} = reshape (msh_side.geo_map(idim,:,:), msh_side.nqn, msh_side.nel); > end > ind = space.sp_patch{iptc}.boundary(iside).dofs; > indices = space.gnum{iptc}(ind); > [grad_eu, F] = sp_eval_msh (u(indices), sp_side, msh_side,'gradient'); > > I=size(F); Lx=zeros(I(2),I(3)); Ly=zeros(I(2),I(3)); > Bx=zeros(I(2),I(3)); By=zeros(I(2),I(3)); > BX=zeros(I(2),I(3)); BY=zeros(I(2),I(3)); > for ix=1:I(2) > for iy=1:I(3) > Lx(ix,iy)=F(1,ix,iy); > Ly(ix,iy)=F(2,ix,iy); > Bx(ix,iy)=eu(2,ix,iy); > By(ix,iy)=eu(1,ix,iy); > end > end > B2=Bx.*Bx+By.*By; > Txx=Bx.*Bx0.5*B2; > FFx=op_F_v (sp_side, msh_side, TFx); %Similar as op_f_v but I only sum all the contribution on a boundary, result is a scalar. > end > FX=FX+FFx; > else > end > end > end > > > > the line > [grad_eu, F] = sp_eval_msh (u(indices), sp_side, msh_side,'gradient'); gives me an error >> Reference to nonexistent field 'shape_function_gradients'. Error in sp_eval_msh (line 153) > space.(field{iopt}) = reshape (space.(field{iopt}), ..., > > > > and I do not understand why, as the same line with 'value' instead of 'gradient' works fine. > Is it a problem of continuity of the gradient on the boundary ? Would it be better to use msh_boundary_from_interior or is it something else. This is because the output structure of sp_eval_boundary_side contains the value of the functions, but not the gradients. I will modify this, to give the possibility of computing also the gradient. In the meantime, it is very easy to modify the code to obtain the gradient. You only need to replace the call to sp_eval_boundary_side by sp_side = sp_precompute (space.sp_patch{iptc}.boundary(iside), msh_side, 'value', true, 'gradient', true) In any case, be aware that sp_side only contains information that can be obtained from the boundary entities. In particular, in the gradient you will only get the _tangential derivatives_. If you want to compute the flux across the boundary, you will need the normal derivative. In this case, you will need to call msh_boundary_side_from_interior. You can see an example of how it is used in the file sp_vector/sp_weak_drchlt_bc(see also Section 4.2.2 in the paper). > Q3. > > In another simulation, I try to solve my original problem in axisymmetric coordinate system with (r,z) instead of (x,y), where y becomes r. > > I had to modify several operators to take into account the additional term in the integrals dxdy becomes rdrdz, so I add input arguments in operators to transport the info problem_data.coord='axi' or 'carth'. > > When building the stiffness matrix and the rhs, I therefore modify the last functions such as : > For the stiffness matrix : in op_gradu_gradv > ... > if strcmp(coord,'axi') > jacdet_weights = msh.jacdet .* msh.quad_weights .* coeff .* y; > else > jacdet_weights = msh.jacdet .* msh.quad_weights .* coeff; > end > > ncounter = 0; > for iel = 1:msh.nel > ... > > where y=x{:,2} is inherited from op_gradu_gradv_tp > > > And the same idea applied for the rhs : in op_f_v > y = reshape (y, spv.ncomp, msh.nqn, msh.nel); > for iel = 1:msh.nel > if (all (msh.jacdet(:,iel))) > if strcmp(coord,'axi') > coeff(:,:,iel) = coeff(:,:,iel).*y(:,:,iel); > end > > ... > > > When computing like this, I still have 30% discrepancy in the simulation result, therefore I tried to add further modifications when imposing the (Dirichlet) boundary conditions and modifying sp_drchlt_l2_proj(both op_u_v_mp and op_f_v_mp). Boundary number 3, the bottom one, is on r = 0, so all the corresponding coefficients in the stiffness matrix M are zero, and the dofs cannot be calculated as the matrix is singular to working precision. > Was it the correct procedure to get to the axisymmetric coordinate? Do you have any further suggestions? I implemented cylindrical coordinates some time ago. The way I did it was slightly different from yours, since I was replacing op_gradu_gradv_tp by a different function (op_gradu_gradv_cyl), where I multiplied by the rcoordinate, instead of doing it inside op_gradu_gradv. But the result should be the same, I don't understand where this 30% of discrepancy is coming from. I suggest you to compute a very simple example with analytical solution in the square, to check where the errors may arise, and to verify that all your data and operators are correctly computed. I hope the answers were helpful to you. Do not hesitate to contact me if you need further help. Best regards, Rafa 
From: Friedrich, L.A.J. <L.Friedrich@tu...>  20170825 13:20:18

Dear GeoPDE users, dear Rafa, I have started using GeoPDE to solve Laplace equations in a 2D multipatch geometry. I have surrounded my geometry with air and imposed null Dirichlet condition on each boundary. I tried to get more familiar with the code in implementing different coefficients for different patches materials, and it worked fine. I am now trying to achieve other tasks, but I encountered several difficulties. I hope you can give me some useful advices and comments on the few questions written below. Q1. In one simulation, I want to impose a source term on an interface between 2 patches, more precisely on top and bottom boundary of one patch. Therefore I used an adapted version of op_f_v_tp where I start to compute the interior rhs by column, the interior points are set to 0. Then I wrote ... for iside = 3 msh_side = msh_eval_boundary_side (msh, iside); sp_side = sp_eval_boundary_side (space, msh_side); for idim = 1:msh.rdim x{idim} = reshape (msh_side.geo_map(idim,:,:), ...msh_side.nqn, msh_side.nel); end coeff_top =f.t(x{:})%inherited from additional input funtions rhs = rhs + op_f_v (sp_side, msh_side, coeff_top); end ... same lines for iside=4, but with coeff_bot as the function. The thing is in one patch I have 104dofs (nsub=[9,4] and nquad=[4,4], degree =[4,4]) which is the same when I compute the rhs by column, but then, on one boundary I have only 13dofs. Therefore how to assemble rhs as they have different dimensions ? Additional question here : we have quadrature points in the interior of the domain, but also on the boundaries right ? is there a way to plot them on top of the geometry? Q2. In another simulation, I want in a postprocessing part, integrate a quantity (deriving from the gradient of my solution) around the boundary of one (group of) patch. Here is the part of code I start to write for one boundary of one patch function [FX] = sp_eval_FX (u,space,msh,geometry,pts,mat) FX=0; for iptc=1:space.npatch if ismember([iptc],mat) for iside = 3 msh_side = msh_eval_boundary_side (msh.msh_patch{iptc}, iside); sp_side = sp_eval_boundary_side (space.sp_patch{iptc}, msh_side); for idim = 1:msh.rdim x{idim} = reshape (msh_side.geo_map(idim,:,:), msh_side.nqn, msh_side.nel); end ind = space.sp_patch{iptc}.boundary(iside).dofs; indices = space.gnum{iptc}(ind); [grad_eu, F] = sp_eval_msh (u(indices), sp_side, msh_side,'gradient'); I=size(F); Lx=zeros(I(2),I(3)); Ly=zeros(I(2),I(3)); Bx=zeros(I(2),I(3)); By=zeros(I(2),I(3)); BX=zeros(I(2),I(3)); BY=zeros(I(2),I(3)); for ix=1:I(2) for iy=1:I(3) Lx(ix,iy)=F(1,ix,iy); Ly(ix,iy)=F(2,ix,iy); Bx(ix,iy)=eu(2,ix,iy); By(ix,iy)=eu(1,ix,iy); end end B2=Bx.*Bx+By.*By; Txx=Bx.*Bx0.5*B2; FFx=op_F_v (sp_side, msh_side, TFx); %Similar as op_f_v but I only sum all the contribution on a boundary, result is a scalar. end FX=FX+FFx; else end end end the line > [grad_eu, F] = sp_eval_msh (u(indices), sp_side, msh_side,'gradient'); gives me an error > Reference to nonexistent field 'shape_function_gradients'. Error in sp_eval_msh (line 153) space.(field{iopt}) = reshape (space.(field{iopt}), ..., and I do not understand why, as the same line with 'value' instead of 'gradient' works fine. Is it a problem of continuity of the gradient on the boundary ? Would it be better to use msh_boundary_from_interior or is it something else. Q3. In another simulation, I try to solve my original problem in axisymmetric coordinate system with (r,z) instead of (x,y), where y becomes r. I had to modify several operators to take into account the additional term in the integrals dxdy becomes rdrdz, so I add input arguments in operators to transport the info problem_data.coord='axi' or 'carth'. When building the stiffness matrix and the rhs, I therefore modify the last functions such as : For the stiffness matrix : in op_gradu_gradv ... if strcmp(coord,'axi') jacdet_weights = msh.jacdet .* msh.quad_weights .* coeff .* y; else jacdet_weights = msh.jacdet .* msh.quad_weights .* coeff; end ncounter = 0; for iel = 1:msh.nel ... where y=x{:,2} is inherited from op_gradu_gradv_tp And the same idea applied for the rhs : in op_f_v y = reshape (y, spv.ncomp, msh.nqn, msh.nel); for iel = 1:msh.nel if (all (msh.jacdet(:,iel))) if strcmp(coord,'axi') coeff(:,:,iel) = coeff(:,:,iel).*y(:,:,iel); end ... When computing like this, I still have 30% discrepancy in the simulation result, therefore I tried to add further modifications when imposing the (Dirichlet) boundary conditions and modifying sp_drchlt_l2_proj(both op_u_v_mp and op_f_v_mp). Boundary number 3, the bottom one, is on r = 0, so all the corresponding coefficients in the stiffness matrix M are zero, and the dofs cannot be calculated as the matrix is singular to working precision. Was it the correct procedure to get to the axisymmetric coordinate? Do you have any further suggestions? I hope those questions were not too long. Thank you for the comments you will give, I have good hope they will increase my comprehension of GeoPDE structure. Best regards, Léo 
From: Rafael Vazquez <vazquez@im...>  20170621 10:29:52

Dear Tony, you are working in the interface space sp_interface = space.sp_patch{patch}.boundary(side) The size of the vector for each interface is rhs_int = zeros(sp_interface.ndof,1); After you compute rhs_int, you have to pass from the local numbering of the interface to the global numbering of the patch (with boundary.dofs), and then from the numbering of the patch to the global multipatch numbering (with gnum) num = space.gnum{patch}(sp_interface.dofs) rhs(num) = rhs(num)  rhs_int; Best regards, Rafa On 21. 06. 17 05:48, tju_zmx wrote: > Dear Rafa > > Now I can get the equivalent load of the interface pressure correctly. > But the problem is that I don't know how to assemble rhs. > How to deal with following code? > > rhs_press = zeros (space.boundary.ndof, 1); > ...... > rhs_press(space.boundary.gnum{iref}) = > rhs_press(space.boundary.gnum{iref}) + op_pn_v (sp_side, msh_side, pval); > ...... > rhs(space.boundary.dofs) = rhs(space.boundary.dofs)  rhs_press; > > Yours Sincerely, > Tony > > 
From: tju_zmx <tju_zmx@12...>  20170621 03:49:14

Dear Rafa Now I can get the equivalent load of the interface pressure correctly. But the problem is that I don't know how to assemble rhs. How to deal with following code? rhs_press = zeros (space.boundary.ndof, 1); ...... rhs_press(space.boundary.gnum{iref}) = rhs_press(space.boundary.gnum{iref}) + op_pn_v (sp_side, msh_side, pval); ...... rhs(space.boundary.dofs) = rhs(space.boundary.dofs)  rhs_press; Yours Sincerely, Tony 
From: Rafael Vazquez <vazquez@im...>  20170620 14:00:59

Dear Tony, there is ongoing research to impose high continuity across patches, either imposing it weakly (with Nitsche's method, for instance), or strongly (see for instance the works by Sangalli and Takacs, Jorg Peters, or Mario Kapl). None of them has been implemented in geopdes, and I am not planning to do it in the near future. About imposing the pressure on the interface, I don't know which kind of problem you want to solve. In any case, you can access the interfaces from msh.msh_patch{patch_number}.boundary(side_number) space.sp_patch{patch_number}.boundary(side_number) Use it at your own risk. ;) Best regards, Rafa On 20. 06. 17 14:15, tju_zmx wrote: > Dear Rafa > > To get the continuous stress, is it possible to control the regularity > of the solution between different patches? > > If I want to apply the pressure on the interface between patches, how > to assemble the righthand side vector correctly? > Do you think it necessary that I create the msh.interface (similar to > msh.boundary)? Is there another way to realise it ? > > Yours Sincerely, > Tony > > 
From: tju_zmx <tju_zmx@12...>  20170620 12:15:52

Dear Rafa To get the continuous stress, is it possible to control the regularity of the solution between different patches? If I want to apply the pressure on the interface between patches, how to assemble the righthand side vector correctly? Do you think it necessary that I create the msh.interface (similar to msh.boundary)? Is there another way to realise it ? Yours Sincerely, Tony 
From: Rafael Vazquez <vazquez@im...>  20170619 15:04:35

Dear Hadi, in this article there is an approximation of the airfoil with Bsplines. http://www.sciencedirect.com/science/article/pii/S0045782514004708?via%3Dihub You will not get an exact representation with Bsplines, but I don't know any work where it has been done with NURBS. Best regards, Rafael 
From: Hadi Minbashian <minbashian@gm...>  20170617 12:54:48

Dear Professor Vazquez, If you could please let me know how I can represent the geometry of an airfoil (e.g. NACA0012) in GeoPDEs either with Bsplines or NURBS. What I have at my hand is a set of points on the boundary of the airfoil and I need to produce a mesh similar to the one depicted in the attached figure in order to do CFD computations. Best regards  Hadi Minbashian (Mr.) PhD Student in Applied Mathematics, Dept. of Mathematics, Faculty of Math. & Comp. Science, Amirkabir University of Technology Hafez Ave., Tehran, Iran Office: Kharazmi Building, Room 338 Tel: (+98)2164542545 Fax: (+98)2166497930 Email: minbashian@... 