At each clue \(n\), we grow all possible \(n\)-ominoes starting from that square. We also record the squares adjacent to the current polyomino. If any of these also contain the number n then the polyomino cannot be part of the solution. Otherwise we introduce a variable to represent that polyomino.
We can then build ZDDs enforcing that each clued square be covered by exactly one polyomino of the appropriate size, ZDDs ensuring every other square is covered by at most one polyomino we enumerated, and lastly ZDDs guaranteeing adjacent polyominoes must differ in size. Taking the intersection of these ZDDs moves us to the next phase.
If the puzzle is sufficiently simple then the solution drops out by now. Otherwise we have partial solutions with some unallocated squares, and the goal is to partition them into polyominoes satisfying the fillomino constraints.
I’m not sure how to proceed with ZDDs. I’d use brute force: for each solution, enumerate all possible ways to carve up the uncovered squares and check if it works. My program simply eliminates solutions where a surrounded uncovered square lies next to a 1, and similarly when a surrounded uncovered 2x1 region lies next to a 2. This seems to completely solve most human-friendly instances, or leave so little unsolved that it is trivial to finish the puzzle by hand.