<commit_before>import numpy as nm

from sfepy.base.base import Struct
from sfepy.terms.terms import Term

class ContactInfo(Struct):
    """
    Various contact-related data of contact terms.
    """
    pass

class ContactTerm(Term):
    r"""
    """
    name = 'dw_contact'
    arg_types = ('material', 'virtual', 'state')
    arg_shapes = {'material' : 'str', 'virtual' : ('D', 'state'), 'state' : 'D'}
    integration = 'surface'

    def __init__(self, *args, **kwargs):
        Term.__init__(self, *args, **kwargs)

        self.ci = None

    @staticmethod
    def function(out, geo, fmode):
        out_qp = nm.zeros((out.shape[0], geo.n_qp) + out.shape[2:],
                          dtype=out.dtype)
        status = geo.integrate(out, nm.ascontiguousarray(out_qp))

        return status

    def get_fargs(self, material, virtual, state,
                  mode=None, term_mode=None, diff_var=None, **kwargs):
        geo, _ = self.get_mapping(virtual)

        region0 = self.region
        region1 = self.region.domain.regions[material]

        if self.ci is None:
            self.ci = ContactInfo()

        print(region0.name)
        print(region0.shape)
        print(region0.facets)
        print(region0.get_facet_indices())
        print(state.field.surface_data[region0.name].fis)

        print(region1.name)
        print(region1.shape)
        print(region1.facets)
        print(region1.get_facet_indices())
        print(state.field.surface_data[region1.name].fis)

        print(geo)
        print(geo.normal)

        mesh_coors = self.region.domain.mesh.coors
        print(mesh_coors[region0.vertices])

        # Uses field connectivity (higher order nodes).
        sd0 = state.field.surface_data[region0.name]
        sd1 = state.field.surface_data[region1.name]

        # Uses mesh connectivity.
        sdg0 = self.region.domain.surface_groups[region0.name]
        sdg1 = self.region.domain.surface_groups[region1.name]

        print(mesh_coors[sdg0.econn])
        print(mesh_coors[sdg1.econn])

        qps = self.get_physical_qps()
        qp_coors = qps.values
        u_qp = self.get(state, 'val').reshape(qp_coors.shape)

        # Deformed QP coordinates.
        coors = u_qp + qp_coors

        print(coors)

        if diff_var is None:
            fmode = 0

        else:
            fmode = 1

        return geo, fmode<commit_msg>update ContactTerm.get_fargs() - test calling ccontres functions - WIP<commit_after>import numpy as nm

from sfepy.base.base import Struct
from sfepy.terms.terms import Term

class ContactInfo(Struct):
    """
    Various contact-related data of contact terms.
    """
    pass

class ContactTerm(Term):
    r"""
    """
    name = 'dw_contact'
    arg_types = ('material', 'virtual', 'state')
    arg_shapes = {'material' : 'str', 'virtual' : ('D', 'state'), 'state' : 'D'}
    integration = 'surface'

    def __init__(self, *args, **kwargs):
        Term.__init__(self, *args, **kwargs)

        self.ci = None

    @staticmethod
    def function(out, geo, fmode):
        out_qp = nm.zeros((out.shape[0], geo.n_qp) + out.shape[2:],
                          dtype=out.dtype)
        status = geo.integrate(out, nm.ascontiguousarray(out_qp))

        return status

    def get_fargs(self, material, virtual, state,
                  mode=None, term_mode=None, diff_var=None, **kwargs):
        geo, _ = self.get_mapping(virtual)

        region0 = self.region
        region1 = self.region.domain.regions[material]

        if self.ci is None:
            self.ci = ContactInfo()

        print(region0.name)
        print(region0.shape)
        print(region0.facets)
        print(region0.get_facet_indices())
        print(state.field.surface_data[region0.name].fis)

        print(region1.name)
        print(region1.shape)
        print(region1.facets)
        print(region1.get_facet_indices())
        print(state.field.surface_data[region1.name].fis)

        print(geo)
        print(geo.normal)

        mesh_coors = self.region.domain.mesh.coors
        print(mesh_coors[region0.vertices])

        # Uses field connectivity (higher order nodes).
        sd0 = state.field.surface_data[region0.name]
        sd1 = state.field.surface_data[region1.name]

        # Uses mesh connectivity.
        sdg0 = self.region.domain.surface_groups[region0.name]
        sdg1 = self.region.domain.surface_groups[region1.name]

        print(mesh_coors[sdg0.econn])
        print(mesh_coors[sdg1.econn])

        qps = self.get_physical_qps()
        qp_coors = qps.values
        u_qp = self.get(state, 'val').reshape(qp_coors.shape)

        # Deformed QP coordinates.
        coors = u_qp + qp_coors

        print(coors)

        ISN = state.field.efaces.T.copy()
        nsd = region0.dim
        ngp = geo.n_qp
        neq = state.n_dof
        nsn = ISN.shape[0]
        nes = ISN.shape[1]
        nen = state.field.gel.n_vertex

        fis0 = region0.get_facet_indices()
        fis1 = region1.get_facet_indices()
        fis = nm.r_[fis0, fis1]
        elementID = fis[:, 0].copy()
        segmentID = fis[:, 1].copy()

        n = len(elementID)
        IEN = state.field.econn
        # Need surface bf, bfg corresponding to field approximation here, not
        # geo...
        H = nm.asfortranarray(geo.bf[0, :, 0, :])
        ps = state.field.gel.surface_facet.poly_space
        gps, gw = self.integral.get_qp(state.field.gel.surface_facet.name)
        bfg = ps.eval_base(gps, diff=True)
        # ?? shape - try in 3D
        dH = nm.asfortranarray(bfg[0, ...])

        X = nm.asfortranarray(state.field.coors)
        Um = nm.asfortranarray(state().reshape((-1, nsd)))
        xx = nm.asfortranarray(X + Um)

        import sfepy.mechanics.extmods.ccontres as cc

        GPs = nm.empty((n*ngp, 2*nsd+6), dtype=nm.float64, order='F')

        # [longestEdge, GPs]=getLongestEdgeAndGPs(n,nsd,ngp,neq,nsn,nes,nen,elementID,segmentID,ISN,IEN,H,xx);
        longestEdge, GPs = cc.get_longest_edge_and_gps(GPs, neq,
                                                       elementID, segmentID,
                                                       ISN, IEN, H, xx)

        print longestEdge
        print X[IEN[elementID[0]]][ISN]

        # gg = state.field.mappings[('Omega', 2, 'volume')]
        # Implement Region.get_diameters(dim).

        #print nm.sqrt(region0.domain.get_element_diameters(state.field.region.cells, gg[0], 0)).max()
        #print (GPs[:, :2] == qp_coors).all()

        #GPs2 = nm.zeros((n*ngp, 2*nsd+6), dtype=nm.float64, order='F')
        #GPs2[:, :2] = qp_coors
        #GPs2[:, 2:4] = nm.repeat(fis, ngp, axis=0)
        #GPs2[:, 4] = nm.finfo(nm.float32).max

        #print (GPs == GPs2).all()
        #print nm.abs(GPs - GPs2).max()

        AABBmin, AABBmax = cc.get_AABB(xx, longestEdge, IEN, ISN,
                                       elementID, segmentID, neq);

        AABBmin = AABBmin - (0.5*longestEdge);
        AABBmax = AABBmax + (0.5*longestEdge);
        N = nm.ceil((AABBmax - AABBmin) / (0.5*longestEdge)).astype(nm.int32)

        print AABBmin, AABBmax, N
        head, next = cc.init_global_search(N, AABBmin, AABBmax, GPs[:,:nsd])

        print head, next

        npd = region0.tdim - 1
        GPs = cc.evaluate_contact_constraints(GPs, ISN, IEN, N,
                                              AABBmin, AABBmax,
                                              head, next, xx,
                                              elementID, segmentID,
                                              npd, neq, longestEdge)
        print GPs[:, 4]

        #
        Gc = nm.zeros(neq, dtype=nm.float64)

        activeGPs = GPs[:, 2*nsd+3]
        epss = 1e9
        keyContactDetection = 1
        keyAssembleKc = 1

        max_num = 4 * nsd * nsn * ngp * GPs.shape[0]
        print 'max_num:', max_num
        vals = nm.empty(max_num, dtype=nm.float64)
        rows = nm.empty(max_num, dtype=nm.int32)
        cols = nm.empty(max_num, dtype=nm.int32)

        aux = cc.assemble_contact_residual_and_stiffness(Gc, vals, rows, cols,
                                                         GPs, ISN, IEN,
                                                         X, Um, H, dH, gw,
                                                         activeGPs, neq, npd,
                                                         epss,
                                                         keyContactDetection,
                                                         keyAssembleKc)
        Gc, vals, rows, cols, num = aux
        print Gc.mean(), num
        print GPs
        print 'true num:', num
        from sfepy.base.base import debug; debug()

        if diff_var is None:
            fmode = 0

        else:
            fmode = 1

        return geo, fmode
{'fraction_non_alphanumeric': 0.11989100817438691, 'fraction_numerical': 0.01407811080835604, 'mean_word_length': 1.8179782469609724, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 25, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '6060285', 'n_tokens_mistral': 3010, 'n_tokens_neox': 2849, 'n_words': 545}