Plotting of Concrete Reinforcement Ratio

About the development of the FEM module/workbench.

Moderator: bernd

User avatar
bernd
Posts: 8486
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: Plotting of Concrete Reinforcement Ratio

Postby bernd » Sat Jul 07, 2018 7:24 am

AFAIK the 13 nodes of the interface face should be returned in both, in the steel box (solid) and in the concrete box (solid).
User avatar
bernd
Posts: 8486
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: Plotting of Concrete Reinforcement Ratio

Postby bernd » Sat Jul 07, 2018 7:25 am

would you post a file to test and to be sure we gone ude the same file for testing
User avatar
bernd
Posts: 8486
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: Plotting of Concrete Reinforcement Ratio

Postby bernd » Sat Jul 07, 2018 7:27 am

deactivate groupmeshing in FEM prefs!
User avatar
HarryvL
Posts: 1052
Joined: Sat Jan 06, 2018 7:38 pm

Re: Plotting of Concrete Reinforcement Ratio

Postby HarryvL » Sat Jul 07, 2018 3:04 pm

bernd wrote:
Sat Jul 07, 2018 7:25 am
would you post a file to test and to be sure we gone ude the same file for testing
Thanks Bernd, here it is:

CON_NONCON.fcstd
(54.39 KiB) Downloaded 17 times

PS: I deactivated group meshing, regenerated GMSH and reran CCX. On return from CCX I print in importToolsFem.py the following for material objects with name "concrete":

Code: Select all

concrete_nodes = femmesh.meshtools.get_femnodes_by_refshape(result_mesh, ref)
print("concrete nodes: {}".format(concrete_nodes))
This prints a list of length 761 ((below in full length :D). Although this may even be truncated, because the print out ends in ", 761L," rather than ", 761L]"

Code: Select all

concrete nodes: [1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L, 23L, 24L, 25L, 26L, 27L, 28L, 29L, 30L, 31L, 32L, 33L, 34L, 35L, 36L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 50L, 51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L, 63L, 64L, 65L, 66L, 67L, 68L, 69L, 70L, 71L, 72L, 73L, 74L, 75L, 76L, 77L, 78L, 79L, 80L, 81L, 82L, 83L, 84L, 85L, 86L, 87L, 88L, 89L, 90L, 91L, 92L, 93L, 94L, 95L, 96L, 97L, 98L, 99L, 100L, 101L, 102L, 103L, 104L, 105L, 106L, 107L, 108L, 109L, 110L, 111L, 112L, 113L, 114L, 115L, 116L, 117L, 118L, 119L, 120L, 193L, 194L, 195L, 196L, 197L, 198L, 199L, 200L, 201L, 202L, 203L, 204L, 205L, 206L, 207L, 208L, 209L, 210L, 211L, 212L, 213L, 214L, 215L, 216L, 217L, 218L, 219L, 220L, 221L, 222L, 223L, 224L, 225L, 226L, 227L, 228L, 229L, 230L, 231L, 232L, 233L, 234L, 235L, 236L, 237L, 238L, 239L, 240L, 241L, 242L, 243L, 244L, 245L, 246L, 247L, 248L, 249L, 250L, 251L, 252L, 253L, 254L, 255L, 256L, 257L, 258L, 259L, 260L, 261L, 262L, 263L, 264L, 265L, 266L, 267L, 268L, 269L, 270L, 271L, 272L, 273L, 274L, 275L, 276L, 277L, 278L, 279L, 280L, 281L, 282L, 283L, 284L, 285L, 286L, 287L, 288L, 289L, 290L, 291L, 292L, 293L, 294L, 295L, 296L, 297L, 298L, 299L, 300L, 301L, 302L, 303L, 304L, 305L, 306L, 307L, 308L, 309L, 310L, 311L, 312L, 313L, 314L, 315L, 316L, 317L, 318L, 319L, 320L, 321L, 322L, 323L, 324L, 325L, 326L, 327L, 328L, 329L, 330L, 331L, 332L, 333L, 334L, 335L, 336L, 337L, 338L, 339L, 340L, 341L, 342L, 343L, 344L, 345L, 346L, 347L, 348L, 349L, 350L, 351L, 352L, 353L, 354L, 355L, 356L, 357L, 358L, 359L, 360L, 361L, 362L, 363L, 364L, 365L, 366L, 367L, 368L, 369L, 370L, 371L, 372L, 373L, 374L, 375L, 376L, 377L, 378L, 379L, 380L, 381L, 382L, 383L, 384L, 385L, 386L, 387L, 388L, 389L, 390L, 391L, 392L, 393L, 394L, 395L, 396L, 397L, 398L, 399L, 400L, 401L, 402L, 403L, 404L, 405L, 406L, 407L, 408L, 409L, 410L, 411L, 412L, 413L, 414L, 415L, 416L, 417L, 418L, 419L, 420L, 421L, 422L, 423L, 424L, 425L, 426L, 427L, 428L, 429L, 430L, 431L, 432L, 433L, 434L, 435L, 436L, 437L, 438L, 439L, 440L, 441L, 442L, 443L, 444L, 445L, 446L, 447L, 448L, 449L, 450L, 451L, 452L, 453L, 454L, 455L, 456L, 457L, 458L, 459L, 460L, 461L, 462L, 463L, 464L, 465L, 466L, 467L, 468L, 469L, 470L, 471L, 472L, 473L, 474L, 475L, 476L, 477L, 478L, 479L, 480L, 481L, 482L, 483L, 484L, 485L, 486L, 487L, 488L, 489L, 490L, 491L, 492L, 493L, 494L, 495L, 496L, 497L, 498L, 499L, 500L, 501L, 502L, 503L, 504L, 505L, 506L, 507L, 508L, 509L, 510L, 511L, 512L, 513L, 514L, 515L, 516L, 517L, 518L, 519L, 520L, 521L, 522L, 523L, 524L, 525L, 526L, 527L, 528L, 529L, 530L, 531L, 532L, 533L, 534L, 535L, 536L, 537L, 538L, 539L, 540L, 541L, 542L, 543L, 544L, 545L, 546L, 547L, 548L, 549L, 550L, 551L, 552L, 553L, 554L, 555L, 556L, 557L, 558L, 559L, 560L, 561L, 562L, 563L, 564L, 565L, 566L, 567L, 568L, 569L, 570L, 571L, 572L, 573L, 574L, 575L, 576L, 577L, 578L, 579L, 580L, 581L, 582L, 583L, 584L, 585L, 586L, 587L, 588L, 589L, 590L, 591L, 592L, 593L, 594L, 595L, 596L, 597L, 598L, 599L, 600L, 601L, 602L, 603L, 604L, 605L, 606L, 607L, 608L, 609L, 610L, 611L, 612L, 613L, 614L, 615L, 616L, 617L, 618L, 619L, 620L, 621L, 622L, 623L, 624L, 625L, 626L, 627L, 628L, 629L, 630L, 631L, 632L, 633L, 634L, 635L, 636L, 637L, 638L, 639L, 640L, 641L, 642L, 643L, 644L, 645L, 646L, 647L, 648L, 649L, 650L, 651L, 652L, 653L, 654L, 655L, 656L, 657L, 658L, 659L, 660L, 661L, 662L, 663L, 664L, 665L, 666L, 667L, 668L, 669L, 670L, 671L, 672L, 673L, 674L, 675L, 676L, 677L, 678L, 679L, 680L, 681L, 682L, 683L, 684L, 685L, 686L, 687L, 688L, 689L, 690L, 691L, 692L, 693L, 694L, 695L, 696L, 697L, 698L, 699L, 700L, 701L, 702L, 703L, 704L, 705L, 706L, 707L, 708L, 709L, 710L, 711L, 712L, 713L, 714L, 715L, 716L, 717L, 718L, 719L, 720L, 721L, 722L, 723L, 724L, 725L, 726L, 727L, 728L, 729L, 730L, 731L, 732L, 733L, 734L, 735L, 736L, 737L, 738L, 739L, 740L, 741L, 742L, 743L, 744L, 745L, 746L, 747L, 748L, 749L, 750L, 751L, 752L, 753L, 754L, 755L, 756L, 757L, 758L, 759L, 760L, 761L,
User avatar
HarryvL
Posts: 1052
Joined: Sat Jan 06, 2018 7:38 pm

Re: Plotting of Concrete Reinforcement Ratio

Postby HarryvL » Sat Jul 07, 2018 3:19 pm

And here is the full code segment for identifying concrete nodes ...
ic[node]==0 means undecided
ic[node]==1 means concrete
ic[node]==2 means not concrete

Code: Select all

                
                for obj in FreeCAD.ActiveDocument.Objects:
                    if obj.isDerivedFrom('App::MaterialObjectPython'):
                        print("object: {}".format(obj))
                        print("object Name: {}".format(obj.Material.get('Name')))
                        print("object.References: {}".format(obj.References))
                        if obj.Material.get('Name') == "Concrete":
                            print("CONCRETE")
                            if obj.References == []:
                                for iic in range(nsr):
                                    if ic[iic] == 0:
                                        ic[iic] = 1
                            else:
                                for ref in obj.References:
                                    concrete_nodes = femmesh.meshtools.get_femnodes_by_refshape(result_mesh, ref)
                                    print(type(concrete_nodes))
                                    print("concrete nodes: {}".format(concrete_nodes))
                                    for cn in concrete_nodes:
                                        ic [cn-1] = 1
                        else:
                            print("NOT CONCRETE")
                            if obj.References == []:
                                for iic in range(nsr):
                                    if ic[iic] == 0:
                                        ic[iic] = 2
                            else:
                                for ref in obj.References:
                                    non_concrete_nodes = femmesh.meshtools.get_femnodes_by_refshape(result_mesh, ref)
                                    for ncn in non_concrete_nodes:
                                        ic [ncn-1] = 2
User avatar
bernd
Posts: 8486
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: Plotting of Concrete Reinforcement Ratio

Postby bernd » Mon Jul 09, 2018 4:27 am

open your file and run the following code ...

Code: Select all

for g in App.ActiveDocument.FEMMeshGmsh.FemMesh.Groups:
    print(App.ActiveDocument.FEMMeshGmsh.FemMesh.getGroupName(g))

I get for your file:

Code: Select all

>>> 
>>> for g in App.ActiveDocument.FEMMeshGmsh.FemMesh.Groups:
...     print(App.ActiveDocument.FEMMeshGmsh.FemMesh.getGroupName(g))
... 
FemConstraintFixed_Nodes
FemConstraintFixed_Faces
FemConstraintPressure_Nodes
FemConstraintPressure_Faces
>>> 
Which seams you have group meshing activated in FEM prefs. If there are groups these groups will be taken. Means to be on the save side in any case deactivate automatic group meshing for analysis reference shapes.

I will have a look at the file anyway ...
User avatar
bernd
Posts: 8486
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: Plotting of Concrete Reinforcement Ratio

Postby bernd » Mon Jul 09, 2018 4:42 am

my file without groups attached ... I opened your file remeshed without groups meshing, deleted all results and rerun ccx.

my (your) code (which only partially runs for me) attached ...

Code: Select all

import femmesh
result_mesh = App.ActiveDocument.Result_mesh.FemMesh

for obj in FreeCAD.ActiveDocument.Objects:
    if obj.isDerivedFrom('App::MaterialObjectPython'):
        print("object: {}".format(obj))
        print("object Name: {}".format(obj.Material.get('Name')))
        print("object.References: {}".format(obj.References))
        if obj.Material.get('Name') == "Concrete":
            print("CONCRETE")
            if obj.References == []:
                for iic in range(nsr):
                    if ic[iic] == 0:
                        ic[iic] = 1
            else:
                for ref in obj.References:
                    concrete_nodes = femmesh.meshtools.get_femnodes_by_refshape(result_mesh, ref)
                    print(type(concrete_nodes))
                    print("concrete nodes: {}".format(concrete_nodes))
                    for cn in concrete_nodes:
                        ic [cn-1] = 1
        else:
            print("NOT CONCRETE")
            if obj.References == []:
                for iic in range(nsr):
                    if ic[iic] == 0:
                        ic[iic] = 2
            else:
                for ref in obj.References:
                    non_concrete_nodes = femmesh.meshtools.get_femnodes_by_refshape(result_mesh, ref)
                    for ncn in non_concrete_nodes:
                        ic [ncn-1] = 2


output:

Code: Select all

object: <App::MaterialObjectPython object>
object Name: Concrete
object.References: [(<Part::PartFeature>, ('Solid1',))]
CONCRETE
  ReferenceShape ... Type: Solid, Object name: BooleanFragments, Object label: BooleanFragments, Element name: Solid1
<type 'list'>
concrete nodes: [1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L, 23L, 24L, 33L, 34L, 35L, 36L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 50L, 51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L, 88L, 89L, 90L, 91L, 92L, 93L, 94L, 95L, 96L, 97L, 98L, 99L, 100L]

Traceback (most recent call last):
  File "<input>", line 18, in <module>
NameError: name 'ic' is not defined
>>> 
how do I get ic, nsr etc ? Are these empty lists I need to create before?


CON_NONCON_bernd1.fcstd
(54.44 KiB) Downloaded 15 times
User avatar
HarryvL
Posts: 1052
Joined: Sat Jan 06, 2018 7:38 pm

Re: Plotting of Concrete Reinforcement Ratio

Postby HarryvL » Mon Jul 09, 2018 9:16 am

bernd wrote:
Mon Jul 09, 2018 4:42 am
... how do I get ic, nsr etc ? Are these empty lists I need to create before?
Thanks for doing this Bernd. The piece of code follows a few lines of existing code defining nsr and -yes- ic also gets initialized, as follows:

Code: Select all

        if 'stress' in result_set:
            stress = result_set['stress']
            nsr=len(stress)
            print("nsr: {}".format(nsr))

            if nsr > 0:
                mstress = []
                prinstress1 = []
                prinstress2 = []
                prinstress3 = []
                shearstress = []
                ps1v = []
                ps2v = []
                ps3v = []
#
#              HarryvL: from here onward new code for handling of concrete reinforcement
#                
                ic=np.zeros(nsr)
#
#               HarryvL: addtional arrays to hold reinforcement ratios and mohr coulomb stress          
#                
                rhx = []
                rhy = []
                rhz = []
                moc = []
#
#               HarryvL: determine concrete / non-concrete nodes
#                
                for obj in FreeCAD.ActiveDocument.Objects:
                    if obj.isDerivedFrom('App::MaterialObjectPython'):
                        print("object: {}".format(obj))
                        print("object Name: {}".format(obj.Material.get('Name')))
                        print("object.References: {}".format(obj.References))
                        if obj.Material.get('Name') == "Concrete":
                            print("CONCRETE")
                            if obj.References == []:
                                for iic in range(nsr):
                                    if ic[iic] == 0:
                                        ic[iic] = 1
                            else:
                                for ref in obj.References:
                                    concrete_nodes = femmesh.meshtools.get_femnodes_by_refshape(result_mesh, ref)
                                    print(type(concrete_nodes))
                                    print("concrete nodes: {}".format(concrete_nodes))
                                    for cn in concrete_nodes:
                                        ic [cn-1] = 1
                        else:
                            print("NOT CONCRETE")
                            if obj.References == []:
                                for iic in range(nsr):
                                    if ic[iic] == 0:
                                        ic[iic] = 2
                            else:
                                for ref in obj.References:
                                    non_concrete_nodes = femmesh.meshtools.get_femnodes_by_refshape(result_mesh, ref)
                                    for ncn in non_concrete_nodes:
                                        ic [ncn-1] = 2
                
                for isv in range(nsr):

                    i=stress.values()[isv]

                    rhox=0.
                    rhoy=0.
                    rhoz=0.
                    mc=0.
                    scxx=i[0]
                    scyy=i[1]
                    sczz=i[2]

                    mstress.append(calculate_von_mises(i))
                    
                    if ic[isv]==1:
#
#                       HarryvL: calculation of reinforcement ratio
#
                        rhox, rhoy, rhoz, scxx, scyy, sczz = calculate_rho(i)

#
#                       HarryvL: for concrete scxx etc. are affected by reinforcement (see calculate_rho(i)). for all other materials scxx etc. are the original stresses
#
                    prin1, prin2, prin3, shear, psv = calculate_principal_stress(i,scxx,scyy,sczz)
                    prinstress1.append(prin1)
                    prinstress2.append(prin2)
                    prinstress3.append(prin3)
                    shearstress.append(shear)
                    ps1v.append(psv[0])
                    ps2v.append(psv[1])
                    ps3v.append(psv[2])

#
#                   reinforcement ratios and mohr coulomb criterion
#
                    rhx.append(rhox)
                    rhy.append(rhoy)
                    rhz.append(rhoz)
                    if ic[isv]==1:
                        mc=calculate_mohr_coulomb(prin1,prin2,prin3)
                    moc.append(mc)
                    
 #                   
 #                from here onward original code for filling result object 
 #                                       
User avatar
HarryvL
Posts: 1052
Joined: Sat Jan 06, 2018 7:38 pm

Re: Plotting of Concrete Reinforcement Ratio

Postby HarryvL » Mon Jul 09, 2018 9:39 am

bernd wrote:
Mon Jul 09, 2018 4:42 am
output:

Code: Select all

object: <App::MaterialObjectPython object>
object Name: Concrete
object.References: [(<Part::PartFeature>, ('Solid1',))]
CONCRETE
  ReferenceShape ... Type: Solid, Object name: BooleanFragments, Object label: BooleanFragments, Element name: Solid1
<type 'list'>
concrete nodes: [1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L, 23L, 24L, 33L, 34L, 35L, 36L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 50L, 51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L, 88L, 89L, 90L, 91L, 92L, 93L, 94L, 95L, 96L, 97L, 98L, 99L, 100L]

Traceback (most recent call last):
  File "<input>", line 18, in <module>
NameError: name 'ic' is not defined
>>> 
Bernd, as I reported here: https://forum.freecadweb.org/viewtopic. ... 80#p242996 I also get 100 nodes. However, that is wrong. It should be (113-13)/2+13=65 nodes. There are 113 nodes, half concrete (Solid1), half steel (Solid2) and 13 shared. The print-out indicates there are 100 nodes in Solid1.
Last edited by HarryvL on Mon Jul 09, 2018 11:01 am, edited 2 times in total.
User avatar
HarryvL
Posts: 1052
Joined: Sat Jan 06, 2018 7:38 pm

Re: Plotting of Concrete Reinforcement Ratio

Postby HarryvL » Mon Jul 09, 2018 9:42 am

HarryvL wrote:
Mon Jul 09, 2018 9:39 am
bernd wrote:
Mon Jul 09, 2018 4:42 am
output:

Code: Select all

object: <App::MaterialObjectPython object>
object Name: Concrete
object.References: [(<Part::PartFeature>, ('Solid1',))]
CONCRETE
  ReferenceShape ... Type: Solid, Object name: BooleanFragments, Object label: BooleanFragments, Element name: Solid1
<type 'list'>
concrete nodes: [1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L, 23L, 24L, 33L, 34L, 35L, 36L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 50L, 51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L, 88L, 89L, 90L, 91L, 92L, 93L, 94L, 95L, 96L, 97L, 98L, 99L, 100L]

Traceback (most recent call last):
  File "<input>", line 18, in <module>
NameError: name 'ic' is not defined
>>> 
Bernd, as I reported here: https://forum.freecadweb.org/viewtopic. ... 80#p242996 I also get 100 nodes after disabling group meshing. However, that is wrong. It should be (113-13)/2+13=65 nodes. There are 113 nodes, half concrete (Solid1), half steel (Solid2) and 13 shared. The print-out indicates there are 100 nodes in Solid1.