~grm08/fluidity/hilbert-for-pyop2

« back to all changes in this revision

Viewing changes to assemble/Geostrophic_Pressure.F90

  • Committer: skramer
  • Date: 2010-11-16 15:37:38 UTC
  • Revision ID: svn-v4:5bf5533e-7014-46e3-b1bb-cce4b9d03719:trunk:2477
Changing the storage of vector_field from vfield%val(dim)%ptr(node) to vfield%val(dim,node).

This was a long desired and discussed change that should improve cache performance and will allow vector fields of more than 3 dimension (internally no vtk support). Just realised it could be done with a few simple seds.

Most important changes need:
- removal of wrap_field() for vector_fields.
- this required some surgery in vtk_write_fields, a horrible piece of code in any case
- traffic and mesh_movement were using wrap_field and now have to do a few more copies.
- halo_update() requires a contiguous buffer in memory, so a copy is now done for halo_update(vfield)
- the wrapped python vector fields available in python_state diagnostics had to be rewired to account for the new storage. This actually simplifies the wrapping a lot.

I realise this is a big commit which might give some trouble when merging in big changes (from for example a branch). The following little script I used might be handy for this. You might only want to apply it to files you have modified. Make sure you back up your entire local checkout before applying any of these:

#/bin/sh
FILES=*/*.F90 */tests/*.F90
# replace any %val(1)%ptr(2) where 2 contains a pair of ()
sed -i 's!%val(\([A-Za-z0-9%_\+\* -]*\))%ptr(\([^)]*([^)]*)[^)]*\))!%val(\1,\2)!g' ${FILES}
# replace any %val(1)%ptr(2) where 2 contains no closing )
sed -i 's!%val(\([A-Za-z0-9%_\+\* -]*\))%ptr(\([^)]*\))!%val(\1,\2)!g' ${FILES}
# replace any %val(1)%ptr not followed by (
sed -i 's!%val(\([A-Za-z0-9%_\+\* -]*\))%ptr!%val(\1,:)!g' ${FILES}
# now svn revert femtools/Sparse_Tools.F90
# as it gives false positives for dcsr_matrices
# luckily it doesn't use vector_fields, so we can just revert


Show diffs side-by-side

added added

removed removed

Lines of Context:
460
460
        ewrite(2, *) "Using " // hpg_name
461
461
        have_hpg = .true.
462
462
        do i = 1, hpg%dim
463
 
          ewrite_minmax(hpg%val(i)%ptr)
 
463
          ewrite_minmax(hpg%val(i,:))
464
464
        end do
465
465
      else
466
466
        ewrite(2, *) "No " // hpg_name
485
485
      assert(ele_count(velocity) == ele_count(gp_rhs))
486
486
    
487
487
      do i = 1, velocity%dim
488
 
        ewrite_minmax(velocity%val(i)%ptr)
 
488
        ewrite_minmax(velocity%val(i,:))
489
489
      end do
490
490
    else
491
491
      velocity => dummy_vector
701
701
    assert(ele_count(positions) == ele_count(mom_rhs))
702
702
 
703
703
    do i = 1, mom_rhs%dim
704
 
      ewrite_minmax(mom_rhs%val(i)%ptr)
 
704
      ewrite_minmax(mom_rhs%val(i,:))
705
705
    end do
706
706
    
707
707
    do i = 1, ele_count(mom_rhs)
711
711
    end do
712
712
    
713
713
    do i = 1, mom_rhs%dim
714
 
      ewrite_minmax(mom_rhs%val(i)%ptr)
 
714
      ewrite_minmax(mom_rhs%val(i,:))
715
715
    end do
716
716
    
717
717
    ewrite(1, *) "Exiting subtract_geostrophic_pressure_gradient"
1587
1587
        call assemble_velocity_ele(i, positions, coriolis, velocity, llump_rhs)
1588
1588
      end do
1589
1589
      do i = 1, coriolis%dim
1590
 
        velocity%val(i)%ptr = velocity%val(i)%ptr / masslump%val
 
1590
        velocity%val(i,:) = velocity%val(i,:) / masslump%val
1591
1591
      end do
1592
1592
    else
1593
1593
      select case(cont)
1606
1606
            call assemble_velocity_ele(i, positions, coriolis, rhs, llump_rhs)
1607
1607
          end do
1608
1608
          do i = 1, rhs%dim
1609
 
            ewrite_minmax(rhs%val(i)%ptr)
 
1609
            ewrite_minmax(rhs%val(i,:))
1610
1610
          end do
1611
1611
          call petsc_solve(velocity, mass, rhs, option_path = solver_path)
1612
1612
          call deallocate(rhs)
1621
1621
    end if
1622
1622
    
1623
1623
    do i = 1, velocity%dim
1624
 
      ewrite_minmax(velocity%val(i)%ptr)
 
1624
      ewrite_minmax(velocity%val(i,:))
1625
1625
    end do
1626
1626
    
1627
1627
    ewrite(1, *) "Exiting velocity_from_coriolis"
1741
1741
        call assemble_coriolis_ele(i, positions, velocity, coriolis, llump_rhs)
1742
1742
      end do
1743
1743
      do i = 1, coriolis%dim
1744
 
        coriolis%val(i)%ptr = coriolis%val(i)%ptr / masslump%val
 
1744
        coriolis%val(i,:) = coriolis%val(i,:) / masslump%val
1745
1745
      end do
1746
1746
    else
1747
1747
      select case(cont)
1772
1772
    end if
1773
1773
    
1774
1774
    do i = 1, coriolis%dim
1775
 
      ewrite_minmax(coriolis%val(i)%ptr)
 
1775
      ewrite_minmax(coriolis%val(i,:))
1776
1776
    end do
1777
1777
    
1778
1778
    ewrite(1, *) "Exiting coriolis_from_velocity"
3009
3009
      & solver_path = trim(base_path) // "/coriolis/coriolis_to_velocity")  
3010
3010
    if(dim == 3) then
3011
3011
      ! Recover the vertical velocity
3012
 
      ewrite_minmax(new_velocity%val(W_)%ptr)
 
3012
      ewrite_minmax(new_velocity%val(W_,:))
3013
3013
      call set(new_velocity, W_, new_w)
3014
 
      ewrite_minmax(new_velocity%val(W_)%ptr)
 
3014
      ewrite_minmax(new_velocity%val(W_,:))
3015
3015
      call deallocate(new_w)
3016
3016
    end if
3017
3017