Here's an example script. I'll assume you use V3.6+
- Create a new color attribute
- Iterate over all materials and
fetch their color
- Iterate over all polygons and fill out the
correct colors into the polygon corners.
Note non-assigned materials or materials with no principled bsdf node will have a black color.
You should be interested in reading bpy.types.Mesh, bpy.types.MeshPolygon, bpy.types.AttributeGroup and bpy.types.Attribute.

import bpy
def main():
obj = bpy.data.objects.get("House")
if obj is None or obj.type != "MESH":
return
mesh = obj.data
color_attribute = mesh.color_attributes.get("Color")
if not color_attribute:
color_attribute = mesh.color_attributes.new("Color", "FLOAT_COLOR", "CORNER")
materials_colors = []
for material_slot in obj.material_slots:
material = material_slot.material
if material is None:
materials_colors.append((0, 0, 0, 0))
continue
principled = next((n for n in material.node_tree.nodes if n.type == "BSDF_PRINCIPLED"), None)
if principled is None:
materials_colors.append((0, 0, 0, 0))
continue
materials_colors.append(principled.inputs["Base Color"].default_value[:])
for polygon in mesh.polygons:
for corner_index in polygon.loop_indices:
color_attribute.data[corner_index].color = materials_colors[polygon.material_index]
if __name__ == "__main__":
main()